MSK SPECTRUM freeze major subset deep dives

library(magrittr)
library(tidyverse)
library(Seurat)
library(readxl)
library(cowplot)
library(colorblindr)
library(viridis)
library(progeny)
library(destiny)

coi <- params$cell_type_super
cell_sort <- params$cell_sort
cell_type_major <- params$cell_type_major
louvain_resolution <- params$louvain_resolution
louvain_cluster <- params$louvain_cluster

1 Cluster markers

1.1 Major T.super markers for cell assign

### load all data ---------------------------------
source("_src/global_vars.R")

# seu_obj <- read_rds(paste0("/work/shah/isabl_data_lake/analyses/16/52/1652/celltypes/", coi, "_processed.rds"))
seu_obj <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/outs_pre/", coi, "_seurat_", louvain_resolution, ".rds"))
# seu_obj <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_highqc.rds"))

myfeatures <- c("umapharmony_1", "umapharmony_2", "sample", louvain_cluster, "doublet", "nCount_RNA", "nFeature_RNA", "percent.mt", "doublet_score")

plot_data_wrapper <- function(cluster_res) {
  cluster_res <- enquo(cluster_res)
  as_tibble(FetchData(seu_obj, myfeatures)) %>% 
    left_join(meta_tbl, by = "sample") %>% 
    rename(cluster = !!cluster_res) %>% 
    mutate(cluster = as.character(cluster),
           tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite)))
}

plot_data <- plot_data_wrapper(louvain_cluster)

1.2 Subtype currated markers

helper_f <- function(x) ifelse(is.na(x), "", x)

markers_v7_super[[coi]] %>% 
  group_by(subtype) %>% 
  mutate(rank = row_number(gene)) %>% 
  spread(subtype, gene) %>% 
  mutate_all(.funs = helper_f) %>% 
  formattable::formattable()

rank	CD4.T.CXCL13	CD4.T.naive	CD4.T.reg	CD4.Th17	CD8.T.CXCL13	CD8.T.cytotoxic	CD8.T.ISG	Cycling.T.NK	dissociated	NK.CD56	NK.cytotoxic
1	CD4	CCR7	CD4	CCR6	CCL3	CCL4	IFI6	ASPM	BTG1	AREG	FCGR3A
2	CD40LG	CD4	FOXP3	IL4I1	CCL4L2	CCL5	IFIT1	CENPF	DNAJB1	FCER1G	FGFBP2
3	CTLA4	CD40LG	IL2RA	IL7R	CD8A	CD8A	IFIT2	HIST1H4C	DUSP1	GNLY	KLRD1
4	CXCL13	IL7R	TNFRSF4	KLRB1	CRTAM	CD8B	IFIT3	HMGB2	DUSP2	KLRC1	KLRF1
5	FKBP5	KLF2	TRAC	LST1	CXCL13	CRTAM	ISG15	MKI67	EGR1	KRT81	PRF1
6	IL6ST	LTB		LTB	FABP5	CST7	MX1	STMN1	EGR3	TRDC	SPON2
7	ITM2A	TCF7		TNFSF13B	GZMB	DTHD1	MX2	TOP2A	FOS	TYROBP
8	MAF	TPT1			HAVCR2	GZMA	RSAD2	TUBA1B	FOSB	XCL1
9	NMB				IFNG	GZMH		TUBB	HSPA1A	XCL2
10	NR3C1				LAG3	GZMK		TYMS	HSPA1B
11	PDCD1				MIR155HG	GZMM			HSPA6
12	TNFRSF4				PHLDA1	HLA-DPB1			JUN
13	TOX2				PTMS	ITM2C			JUNB
14	TSHZ2				RBPJ	KLRG1			JUND
15					TNFRSF9	TRGC2			KLF2
16									KLF6
17									MT1E
18									MT1X
19									NR4A1
20									ZFP36L1

1.3 Subtype cluster markers

## define patient specific clusters
patient_specific_clusters <- plot_data %>% 
  group_by(patient_id_short, cluster) %>% 
  tally %>% 
  group_by(cluster) %>% 
  mutate(nrel = n / sum(n),
         ntotal = sum(n)) %>% 
  arrange(desc(nrel)) %>% 
  distinct(cluster, .keep_all = T) %>% 
  filter(nrel > 0.5) %>% 
  ungroup %>% 
  mutate(cluster_label_ps = make.names(paste0("OV.", patient_id_short, ".specific"), unique = T),
         cluster = as.numeric(cluster)) %>% 
  distinct(cluster, cluster_label_ps)



## Hypergeometric test --------------------------------------

# marker_tbl <- read_tsv(paste0("/work/shah/isabl_data_lake/analyses/16/52/1652/celltypes/", coi, "_markers.tsv")) %>% 
#   filter(resolution == louvain_resolution)
marker_tbl <- read_tsv(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/outs_pre/", coi, "_markers_", louvain_resolution, ".tsv"))
# marker_tbl <- read_tsv(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_highqc_markers_02.tsv"))

test_set <- marker_tbl %>% 
  group_by(cluster) %>% 
  filter(!str_detect(gene, "^RPS|^RPL")) %>% 
  slice(1:50) %>% 
  mutate(k = length(cluster)) %>% 
  ungroup %>%
  select(cluster, gene, k) %>% 
  mutate(join_helper = 1) %>% 
  group_by(cluster, join_helper, k) %>%
  nest(test_set = gene)

markers_doub_tbl <- markers_v7 %>% 
  enframe("subtype", "gene") %>% 
  filter(!(subtype %in% unique(c(coi, cell_type_major)))) %>% 
  unnest(gene) %>% 
  group_by(gene) %>% 
  filter(length(gene) == 1) %>% 
  mutate(subtype = paste0("doublet.", subtype)) %>% 
  bind_rows(tibble(subtype = "Mito.high", gene = grep("^MT-", rownames(seu_obj), value = T)))

ref_set <- markers_v7_super[[coi]] %>% 
  bind_rows(markers_doub_tbl) %>% 
  group_by(subtype) %>% 
  mutate(m = length(gene),
         n = length(rownames(seu_obj))-m,
         join_helper = 1) %>% 
  group_by(subtype, m, n, join_helper) %>%
  nest(ref_set = gene)

hyper_tbl <- test_set %>% 
  left_join(ref_set, by = "join_helper") %>% 
  group_by(cluster, subtype, m, n, k) %>%
  do(q = length(intersect(unlist(.$ref_set), unlist(.$test_set)))) %>%
  mutate(pval = 1-phyper(q = q, m = m, n = n, k = k)) %>%
  ungroup %>%
  mutate(qval = p.adjust(pval, "BH"),
         sig = qval < 0.01)

# hyper_tbl %>% 
#   group_by(subtype) %>% 
#   filter(any(qval < 0.01)) %>%
#   ggplot(aes(subtype, -log10(qval), fill = sig)) +
#   geom_bar(stat = "identity") +
#   facet_wrap(~cluster) +
#   coord_flip()
  
low_rank <- str_detect(unique(hyper_tbl$subtype), "doublet|dissociated")
subtype_lvl <- c(sort(unique(hyper_tbl$subtype)[!low_rank]), sort(unique(hyper_tbl$subtype)[low_rank]))
  
cluster_label_tbl <- hyper_tbl %>% 
  mutate(subtype = ordered(subtype, levels = subtype_lvl)) %>% 
  arrange(qval, subtype) %>%
  group_by(cluster) %>% 
  slice(1) %>% 
  mutate(subtype = ifelse(sig, as.character(subtype), paste0("unknown_", cluster))) %>% 
  select(cluster, cluster_label = subtype) %>% 
  ungroup %>% 
  mutate(cluster_label = make.unique(cluster_label, sep = "_")) %>% 
  left_join(patient_specific_clusters, by = "cluster") %>% 
  mutate(cluster_label = ifelse(is.na(cluster_label_ps), cluster_label, cluster_label_ps)) %>% 
  select(-cluster_label_ps)

seu_obj$cluster_label <- unname(deframe(cluster_label_tbl)[as.character(unlist(seu_obj[[paste0("RNA_snn_res.", louvain_resolution)]]))])
plot_data$cluster_label <- seu_obj$cluster_label

cluster_n_tbl <- seu_obj$cluster_label %>% 
  table() %>% 
  enframe("cluster_label", "cluster_n") %>% 
  mutate(cluster_nrel = cluster_n/sum(cluster_n))

marker_sheet <- marker_tbl %>% 
  left_join(cluster_label_tbl, by = "cluster") %>% 
  group_by(cluster_label) %>% 
  filter(!str_detect(gene, "^RPS|^RPL")) %>% 
  slice(1:50) %>% 
  mutate(rank = row_number(-avg_logFC)) %>% 
  select(cluster_label, gene, rank) %>% 
  ungroup %>% 
  mutate(cluster_label = ordered(cluster_label, levels = unique(c(names(clrs$cluster_label[[coi]]), sort(cluster_label))))) %>% 
  spread(cluster_label, gene) %>% 
  mutate_all(.funs = helper_f)

marker_tbl_annotated <- marker_tbl %>% 
  left_join(cluster_label_tbl, by = "cluster") %>% 
  left_join(cluster_n_tbl, by = "cluster_label") %>% 
  select(-cluster, -resolution) %>% 
  mutate(cluster_label = ordered(cluster_label, levels = unique(c(names(clrs$cluster_label[[coi]]), sort(cluster_label))))) %>% 
  arrange(cluster_label, -avg_logFC, p_val_adj)

write_tsv(marker_sheet, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_marker_sheet.tsv"))

write_tsv(marker_sheet, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/supplementary_tables/", coi, "_marker_sheet.tsv"))

write_tsv(marker_tbl_annotated, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_marker_table_annotated.tsv"))

write_tsv(marker_tbl_annotated, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/supplementary_tables/", coi, "_marker_table_annotated.tsv"))

formattable::formattable(marker_sheet)

rank	CD4.T.naive	CD4.T.CXCL13	CD4.T.reg	CD4.Th17	CD8.T.cytotoxic	CD8.T.CXCL13	CD8.T.ISG	NK.CD56	NK.cytotoxic	Cycling.T.NK	dissociated	dissociated_1	doublet.Monocyte	doublet.Ovarian.cancer.cell	Mito.high	OV.008.specific	OV.008.specific.1
1	IL7R	CXCL13	TNFRSF4	KLRB1	GZMK	CXCL13	IFIT3	GNLY	FGFBP2	STMN1	CCL4L2	HSPA1A	HLA-DRA	WFDC2	XIST	KLF2	SOX4
2	CCR7	NMB	IL2RA	IL7R	CD8A	GZMB	ISG15	TYROBP	FCGR3A	MKI67	CCL4	HSPA1B	CST3	SLPI	VMP1	CCR7	PTCRA
3	KLF2	NR3C1	FOXP3	IL4I1	CD8B	CCL4L2	IFIT1	AREG	SPON2	TUBA1B	IFNG	MT1X	SPP1	KRT18	NEAT1	JUNB	MAL
4	EEF1B2	FKBP5	CTLA4	CEBPD	ITM2C	MIR155HG	MX1	KLRC1	PRF1	TOP2A	FOS	DNAJB1	CXCL8	LCN2	MALAT1	SELL	MZB1
5	TPT1	IL6ST	LTB	LST1	GZMH	TNFRSF9	IFIT2	FCER1G	KLRF1	CENPF	FOSB	MT1E	S100A9	KRT19	IL7R	DUSP1	DNTT
6	EEF1A1	MAF	RTKN2	TNFSF13B	CCL5	HAVCR2	RSAD2	TRDC	GNLY	HIST1H4C	TNF	HSPA6	S100A8	S100A13	SPOCK2	FOS	TFDP2
7	TCF7	ITM2A	BATF	LTB	TRGC2	RBPJ	IFI6	KRT81	KLRD1	TUBB	CCL3	HSP90AA1	CD74	MDK	MTRNR2L12	IL7R	STMN1
8	MAL	CTLA4	TNFRSF18	FOS	GZMA	LAG3	MX2	XCL1	NKG7	HMGB2	CCL3L1	HSPH1	HLA-DQA1	KRT8	NABP1	EEF1A1	AC084033.3
9	CD40LG	TSHZ2	SAT1	SLC4A10	KLRG1	IFNG	ISG20	KLRD1	CX3CR1	ASPM	EGR1	HSPE1	FTL	CLU	PBXIP1	BTG2	CD1E
10	GPR183	LIMS1	TBC1D4	CCR6	CCL4	PTMS	IFI44L	IGFBP2	GZMB	TYMS	JUN	HSPD1	APOE	SPINT2	AC058791.1	CD69	ARPP21
11	LDHB	CD40LG	TIGIT	NFKBIA	CRTAM	CCL3	HERC5	XCL2	PLAC8	NUSAP1	NFKBID	HSPB1	HLA-DQB1	CALD1	CD4	GPR183	CDK6
12	SELL	PDCD1	GADD45A	RORA	CST7	CD8A	OAS1	CLIC3	CLIC3	HMGN2	NR4A2	HSP90AB1	LYZ	CD24	MIAT	PIK3IP1	MAP1A
13	SNHG8	TNFRSF4	TNFRSF1B	NCR3	GZMM	CRTAM	SAMD9L	KRT86	PLEK	PCLAF	EGR2	DNAJA1	HLA-DPB1	IFITM3	LTB	DNAJB1	GLUL
14	NOSIP	CD4	PMAIP1	TNFAIP3	DTHD1	FABP5	TNFSF10	IL2RB	TYROBP	H2AFZ	CD69	CACYBP	MARCKS	VEGFA	SELL	CD55	AC011893.1
15	PABPC1	RNF19A	UGP2	AQP3	HLA-DPB1	PHLDA1	STAT1	CEBPD	PTGDS	SMC4	XCL2	JUN	AIF1	SOX4	MT-ND3	RACK1	GRASP
16	NOP53	RBPJ	IKZF2	MYBL1	PPP1R14B	JAML	EIF2AK2	CTSW	EFHD2	PCNA	TNFSF9	HSPA8	C1QB	C19orf33	CCR7	TSC22D3	ADA
17	LEF1	CORO1B	TNFRSF9	DUSP1	CD3G	TIGIT	GBP1	TXK	FCER1G	HIST1H1B	IER2	FKBP4	C1QA	TM4SF1	MT-ND6	LDHB	AL357060.1
18	LTB	CPM	ICOS	LINC01871	THEMIS	KRT86	MT2A	MATK	CST7	TPX2	AC020916.1	CHORDC1	FN1	RAB13	MT-ND1	FKBP5	CD1B
19	EIF3E	ZBED2	SOX4	CCL20	EOMES	CCL5	OAS3	KLRB1	GZMH	UBE2C	NR4A1	RGS2	FTH1	MAL2	RORA	SARAF	MIR181A1HG
20	RACK1	AHI1	LINC01943	JAML	TC2N	CXCR6	EPSTI1	CD7	ADGRG1	DUT	DUSP2	DNAJB4	BASP1	YBX3	CD28	EEF1B2	CCDC26
21	JUNB	ICA1	IL32	CTSH	DUSP2	LINC01871	SAMD9	CD63	CCL3	CLSPN	KLF6	DNAJA4	C15orf48	ELF3	MT-CO3	TCF7	CLDN5
22	NACA	TOX2	ARID5B	TMIGD2	LYAR	PDCD1	PLSCR1	NKG7	HOPX	SMC2	GADD45B	ZFAND2A	APOC1	TNFAIP2	MT-CYB	CXCR4	SOCS2
23	SOCS3	DUSP4	LAYN	TNFRSF25	CD3D	HLA-DRB1	IFI35	CCL3	ZEB2	CKS1B	DUSP1	ANXA1	G0S2	COL1A1	MT-ATP6	PLAC8	RCAN1
24	TOMM7	AC004585.1	CD27	ERN1	PPP2R5C	TNIP3	XAF1	TMIGD2	IGFBP7	UBE2S	BTG2	PPP1R15A	MNDA	MT1G	MT-CO1	TPT1	JCHAIN
25	UBA52	ARID5B	BIRC3	DPP4	SLF1	GAPDH	OASL	HOPX	PRSS23	ATAD2	CRTAM	DUSP1	NPC2	KRT7	N4BP2L2	AP3M2	ID1
26	TMEM123	CCDC50	CORO1B	KLRG1	KIAA1551	CD63	IFI44	TNFRSF18	CD247	TMPO	ZFP36L1	FOS	MS4A6A	CST3	PTPRC	PPP1R15A	VIPR2
27	SERINC5	CD84	TYMP	SPOCK2	CXCR6	FAM3C	USP18	CMC1	AKR1C3	TUBB4B	PPP1R15A	FOSB	LST1	IFI27	SF1	CCND3	CASC15
28	EEF2	IGFL2	CD4	ZBTB16	YBX3	GZMH	LY6E	GSTP1	C1orf21	PTTG1	FASLG	SERPINH1	IL1B	CLDN3	PCSK7	NACA	GALNT2
29	FXYD5	SRGN	DUSP4	IFNGR1	HLA-DPA1	CTLA4	CMPK2	SRGAP3	MYBL1	KNL1	TAGAP	UBC	GSN	BCAM	MT-ND2	NOP53	HES4
30	TRABD2A	BATF	ENTPD1	FKBP11	CCL4L2	CCND2	NT5C3A	KLRC2	AREG	CENPE	NR4A3	AHSA1	GRN	CCND1	MT-CO2	ZBTB16	PFKFB2
31	TSHZ2	CH25H	CTSC	TPT1	F2R	SPRY1	DDX58	LAT2	S1PR5	RRM2	RGCC	TSPYL2	CD83	S100A1	AAK1	SC5D	MARCKSL1
32	SARAF	TNFRSF18	MIR4435-2HG	MGAT4A	CXCR4	CD8B	HELZ2	GZMB	CTSW	DEK	KDM6B	UBB	CTSH	RBP1	TSHZ2	NOSIP	AC002454.1
33	AQP3	SPOCK2	LINC02099	PDCD4	FAM102A	VCAM1	IRF7	LINC00996	KLF2	HMGB1	JUND	NEU1	MEF2C	NUPR1	NKTR	FOSB	APBA2
34	ANK3	CHN1	MAGEH1	S100A4	HLA-DRB1	HLA-DRA	TRIM22	NCAM1	ABHD17A	BIRC5	TNFAIP3	DNAJB6	CYBB	S100A16	IKZF1	ERAP2	NREP
35	RIPOR2	CD200	SPOCK2	SATB1	SLAMF7	ID2	PARP14	PRF1	TTC38	HELLS	EGR3	GADD45B	CTSB	SNCG	PIK3IP1	PLK3	TSHR
36	AP3M2	RGS1	CARD16	S100A6	APOBEC3G	ITGAE	DDX60	CXXC5	PTPN12	DLGAP5	TNFSF14	BTG2	C1QC	DSP	MT-ND5	EEF1D	TP53INP1
37	TOB1	RILPL2	PHACTR2	CD40LG	SH2D1A	DUSP4	LAG3	IFITM3	KLRB1	CKS2	SLA	BAG3	CD14	TCEAL9	CD2	BTG1	YBX3
38	ZFAS1	ZNRF1	S100A4	B3GALT2	CD84	LYST	DDX60L	ZNF683	CCL4	NASP	DUSP6	KLF6	GLUL	SELENOP	IL6R	ZFP36L2	MME
39	LINC02273	TNFRSF25	STAM	ABCB1	STK17A	ENTPD1	IFIH1	KLRF1	PTGDR	H2AFX	GZMK	ERN1	CSF3R	EPCAM	B2M	EEF2	CHI3L2
40	EIF4B	METTL8	SPATS2L	RUNX2	CCR5	SRGAP3	PPM1K	MCTP2	TRDC	CDK1	DDX3X	JUNB	SGK1	MUC1	CCR4	HNRNPA1	LEF1
41	ANXA1	SLA	GLRX	TLE1	TUBA4A	TNFSF4	OAS2	SH2D1B	ITGB2	H2AFV	ZFP36	JUND	SPI1	CTTN	RNF213	FBXO32	UHRF1
42	NSA2	SMCO4	AC005224.3	EEF1A1	ARAP2	NDFIP2	RNF213	ITGA1	XBP1	FABP5	RASGEF1B	DEDD2	FCGR2A	ITGB8	SYNE2	PABPC1	SSBP2
43	TNFRSF25	BTLA	MAF	CERK	GPR174	GOLIM4	PARP9	IFITM2	CEP78	CDKN3	REL	IER5L	MS4A7	CRABP2	MT-ND4	PASK	ARRDC2
44	CTSL	SESN3	BTG3	RORC	PECAM1	AKAP5	PNPT1	CCL5	CMC1	CXCL13	KLRG1	CLK1	BCL2A1	C3	DDX17	VSIR	BCL11A
45	SESN3	NAP1L4	PBXIP1	PERP		CD27	SP110	CD38	LITAF	GAPDH	TSC22D3	H3F3B	PLAUR	S100A14	PNISR	EIF3H	SMIM3
46	EEF1D	BHLHE40-AS1	F5	LTK		HLA-DPA1	SAT1	SLC16A3	BIN2	MCM7	CD8A	NR4A1	RNASE1	CKB	STK4	TXNIP	LRRC28
47	FAU	MIR155HG	SLAMF1	PLCB1		SNAP47	C19orf66	ITGAX	CHST2	EZH2	MCL1	CD69	MARCH1	CLDN7	CD44	LDLRAP1	GNA15
48	LDLRAP1	BIRC3	IL1R1	LTC4S		RGS1	STAT2	CAPN12	CD300A	GTSE1	FAM53C	DDIT3	CPVL	IGF2BP2	CD3D	CMTM8	RUFY3
49	FLT3LG	PTPN13	DNPH1	EEF1B2		ITM2A	BST2	CD247	ARL4C	ZWINT	CCL5	EIF4A2	SERPINA1	SOX17	HLA-A	SCML1	SCAI
50	TIMP1	CYSLTR1	TRAC	KIF5C		APOBEC3C	LAP3	SAMD3	TXK	RANBP1	CXCR4	CXCR4	ALDH2	SMIM22	HLA-E	LINC00402	ATP6AP1L

2 Clusters

2.1 sizes

enframe(sort(table(seu_obj$cluster_label))) %>% 
  mutate(name = ordered(name, levels = rev(name))) %>% 
  ggplot() +
  geom_bar(aes(name, value), stat = "identity") +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
  labs(y = c("#cells"), x = "cluster")

2.2 UMAP

alpha_lvl <- ifelse(nrow(plot_data) < 20000, 0.2, 0.1)
pt_size <- ifelse(nrow(plot_data) < 20000, 0.2, 0.05)

common_layers_disc <- list(  
  geom_point(size = pt_size, alpha = alpha_lvl),
  NoAxes(),
  guides(color = guide_legend(override.aes = list(size = 2, alpha = 1))),
  labs(color = "")
)

common_layers_cont <- list(  
  geom_point(size = pt_size, alpha = alpha_lvl),
  NoAxes(),
  scale_color_gradientn(colors = viridis(9)),
  guides(color = guide_colorbar())
)

ggplot(plot_data, aes(umapharmony_1, umapharmony_2, color = cluster_label)) + 
  common_layers_disc +
  #facet_wrap(~therapy) +
  ggtitle("Sub cluster")

3 Filtering out doublet clusters

3.1 UMAP

my_subtypes <- names(clrs$cluster_label[[coi]])
my_subtypes <- c(my_subtypes, unlist(lapply(paste0("_", 1:3), function(x) paste0(my_subtypes, x)))) %>% .[!str_detect(., "doublet|dissociated|Mito.high")]
my_subtypes <- my_subtypes[my_subtypes %in% unique(seu_obj$cluster_label)]
my_subtypes <- my_subtypes[my_subtypes %in% names(clrs$cluster_label[[coi]])]

if (cell_sort == "CD45+") {
cells_to_keep <- colnames(seu_obj)[seu_obj$cluster_label %in% my_subtypes & !(str_detect(seu_obj$cell_id, "CD45N"))]
}

if (cell_sort == "CD45-") {
cells_to_keep <- colnames(seu_obj)[seu_obj$cluster_label %in% my_subtypes & !(str_detect(seu_obj$cell_id, "CD45P"))]
}

# seu_obj_sub <- subset(seu_obj, cells = cells_to_keep)
# seu_obj_sub <- RunUMAP(seu_obj_sub, dims = 1:50, reduction = "harmony", reduction.name = "umapharmony", reduction.key = "umapharmony_")
# seu_obj_sub$cluster_label <- seu_obj$cluster_label[colnames(seu_obj) %in% colnames(seu_obj_sub)]
# write_rds(seu_obj_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered.rds"))
seu_obj_sub <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered.rds"))

plot_data_sub <- as_tibble(FetchData(seu_obj_sub, c(myfeatures, "cluster_label"))) %>% 
  left_join(meta_tbl, by = "sample") %>% 
  mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
         tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite))) %>% 
  mutate(cell_id = colnames(seu_obj_sub),
         cluster_label = ordered(cluster_label, levels = my_subtypes),
         )
  
ggplot(plot_data_sub, aes(umapharmony_1, umapharmony_2, color = cluster_label)) + 
  common_layers_disc +
  ggtitle("Sub cluster") +
  #facet_wrap(~cluster_label) +
  scale_color_manual(values = clrs$cluster_label[[coi]])

ggplot(plot_data_sub, aes(umapharmony_1, umapharmony_2, color = patient_id_short)) + 
  common_layers_disc +
  ggtitle("Patient") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$patient_id_short)

ggplot(plot_data_sub, aes(umapharmony_1, umapharmony_2, color = tumor_supersite)) + 
  # geom_point(aes(umapharmony_1, umapharmony_2), 
  #            color = "grey90", size = 0.01, 
  #            data = select(plot_data_sub, -tumor_supersite)) +
  common_layers_disc +
  ggtitle("Site") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$tumor_supersite)

write_tsv(select(plot_data_sub, cell_id, everything(), -umapharmony_1, -umapharmony_2, -contains("RNA_")), paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_embedding.tsv"))

4 sub clustering

big_clusters <- list(
  CD4.T.CXCL13 = "CD4.T.CXCL13", 
  CD4.T.naive = "CD4.T.naive", 
  CD4.T.reg = "CD4.T.reg", 
  CD4.Th17 = "CD4.Th17", 
  CD8.T.CXCL13 = "CD8.T.CXCL13", 
  CD8.T.cytotoxic = "CD8.T.cytotoxic", 
  CD8.T.ISG = "CD8.T.ISG", 
  Cycling.T.NK = "Cycling.T.NK", 
  NK.CD56 = "NK.CD56", 
  NK.cytotoxic = "NK.cytotoxic"
)

cells_to_keep_list <- lapply(big_clusters, function(x) 
colnames(seu_obj_sub)[seu_obj_sub$cluster_label %in% x])

# seu_list <- lapply(cells_to_keep_list, function(x) subset(seu_obj_sub, cells = x))
# 
# preprocess_wrapper <- . %>%
#   FindNeighbors(reduction = "harmony", dims = 1:50) %>%
#   # FindClusters(res = 0.1) %>%
#   # FindClusters(res = 0.2) %>%
#   FindClusters(res = 0.3)
# 
# seu_list <- lapply(seu_list, preprocess_wrapper)
# write_rds(seu_list, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered_sub_clusters.rds"))

seu_list <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered_sub_clusters.rds"))

get_cluster_vector <- function(cluster_label, cluster_res = RNA_snn_res.0.3) {
  cluster_res <- enquo(cluster_res)
  cbind(cell_id = colnames(seu_list[[cluster_label]]), FetchData(seu_list[[cluster_label]], c(as_label(cluster_res)))) %>%
    as_tibble %>%
    mutate(cluster_extended = paste0(cluster_label, "_", !!cluster_res)) %>%
    select(cell_id, cluster_extended) %>% 
    deframe
}

cluster_vector_list <- lapply(big_clusters, get_cluster_vector)

cluster_vector <- unlist(cluster_vector_list, use.names = F) %>% setNames(lapply(cluster_vector_list, names) %>% unlist(use.names = F))

seu_obj$cluster_extended <- cluster_vector[seu_obj$cell_id] %>%
  setNames(seu_obj$cell_id)

seu_obj_sub$cluster_extended <- cluster_vector[seu_obj_sub$cell_id] %>%
  setNames(seu_obj_sub$cell_id)

cluster_extended_uniq <- sort(as.character(na.omit(unique(seu_obj_sub$cluster_extended))))

Idents(seu_obj_sub) <- seu_obj_sub$cluster_extended
Idents(seu_obj) <- seu_obj$cluster_extended

patient_specific_subcluster <- FetchData(seu_obj_sub, c("patient_id", "cluster_extended")) %>% 
  as_tibble() %>% 
  group_by(patient_id, cluster_extended) %>% 
  tally %>% 
  group_by(cluster_extended) %>% 
  mutate(nrel = n/sum(n),
         ntotal = sum(n)) %>% 
  arrange(desc(nrel)) %>% 
  distinct(cluster_extended, .keep_all = T) %>% 
  ungroup

# marker_tbl_extended <- lapply(
#   cluster_extended_uniq,
#   function(x) {
#     marker_tbl <- FindMarkers(seu_obj_sub, ident.1 = x) %>%
#       as_tibble(rownames = "gene") %>%
#       separate(gene, into = c("gene", "drop"), sep = "\\.\\.\\.") %>%
#       select(-drop) %>%
#       mutate(cluster_extended = x)
#     write_tsv(marker_tbl, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_marker_table_subcluster_", x, ".tsv"))
#     return(marker_tbl)
# }) %>% bind_rows
# 
# write_tsv(marker_tbl_extended, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_marker_table_annotated_extended.tsv"))

marker_tbl_extended <- read_tsv(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_marker_table_annotated_extended.tsv"))

cluster_label_extended <- c(
  CD4.T.naive_0 = "CD4.T.effector.memory", 
  CD4.T.naive_1 = "CD4.T.naive.centr.mem.2", ############## updated in v7
  CD4.T.naive_2 = "CD4.T.naive.centr.mem.1", ############## updated in v7
  CD4.T.CXCL13_0 = "CD4.T.dysfunc.early", 
  CD4.T.CXCL13_1 = "CD4.T.dysfunc.late.1", 
  CD4.T.CXCL13_2 = "CD4.T.dysfunc.late.2", 
  CD4.T.reg_0 = "CD4.T.reg.2", 
  CD4.T.reg_1 = "CD4.T.reg.1", 
  CD4.T.reg_2 = "CD4.T.reg.ISG", 
  CD4.T.reg_3 = "CD4.T.reg.3", 
  CD4.Th17_0 = "CD4.Th17.1", 
  CD4.Th17_1 = "CD4.Th17.2", 
  CD4.Th17_2 = "CD4.Th17.3", ############## updated in v7
  CD8.T.cytotoxic_0 = "CD8.T.effector.memory", 
  CD8.T.cytotoxic_1 = "CD8.T.cytotoxic", 
  CD8.T.cytotoxic_2 = "dissociated_3", 
  CD8.T.cytotoxic_3 = "CD8.T.naive.centr.mem", 
  CD8.T.cytotoxic_4 = "gd.T.cell", 
  CD8.T.cytotoxic_5 = "doublet.Plasma.cell_1", 
  CD8.T.CXCL13_0 = "CD8.T.dysfunc.early", 
  CD8.T.CXCL13_1 = "CD8.T.dysfunc.late", 
  CD8.T.CXCL13_2 = "NK.reg.CRTAM", 
  CD8.T.ISG_0 = "CD8.T.ISG.early", ############# updated in v7
  CD8.T.ISG_1 = "CD4.T.ISG", 
  CD8.T.ISG_2 = "CD8.T.dysfunc.ISG", ############# updated in v7
  CD8.T.ISG_3 = "CD8.T.ISG.late", ############# updated in v7
  CD8.T.ISG_4 = "OV.009.specific", ############## updated in v7
  NK.CD56_0 = "NK.reg.IGFBP2",
  NK.CD56_1 = "NK.reg.CD56", 
  NK.CD56_2 = "NK.reg.ISG", 
  NK.CD56_3 = "NK.reg.KRT81.KRT86.1", ############## updated in v7
  NK.CD56_4 = "NK.reg.KRT81.KRT86.2", ############## updated in v7
  NK.CD56_5 = "NK.reg.CCL3",
  NK.cytotoxic_0 = "NK.cytotoxic.SPON2.2", 
  NK.cytotoxic_1 = "NK.cytotoxic.GZMH", 
  NK.cytotoxic_2 = "NK.cytotoxic.SPON2.1",
  NK.cytotoxic_3 = "dissociated_4", ############## updated in v7
  Cycling.T.NK_0 = "Cycling.CD8.T.1",
  Cycling.T.NK_1 = "Cycling.CD8.T.2", 
  Cycling.T.NK_2 = "Cycling.CD8.T.3", 
  Cycling.T.NK_3 = "Cycling.CD8.T.4", 
  Cycling.T.NK_4 = "Cycling.NK.3", 
  Cycling.T.NK_5 = "Cycling.CD4.T", 
  Cycling.T.NK_6 = "Cycling.NK.2", 
  Cycling.T.NK_7 = "Cycling.NK.1"
)

names(seu_obj$cluster_label) <- colnames(seu_obj)
seu_obj$cluster_label[is.na(seu_obj$cluster_label)] <- "NA"
Idents(seu_obj) <- seu_obj$cluster_label
seu_obj$cluster_label_sub <- cluster_label_extended[seu_obj$cluster_extended]
names(seu_obj$cluster_label_sub) <- colnames(seu_obj)
seu_obj$cluster_label_sub[is.na(seu_obj$cluster_label_sub)] <- "NA"
Idents(seu_obj) <- seu_obj$cluster_label_sub

names(seu_obj_sub$cluster_label) <- colnames(seu_obj_sub)
seu_obj_sub$cluster_label[is.na(seu_obj_sub$cluster_label)] <- "NA"
Idents(seu_obj_sub) <- seu_obj_sub$cluster_label
seu_obj_sub$cluster_label_sub <- cluster_label_extended[seu_obj_sub$cluster_extended]
names(seu_obj_sub$cluster_label_sub) <- colnames(seu_obj_sub)
seu_obj_sub$cluster_label_sub[is.na(seu_obj_sub$cluster_label_sub)] <- "NA"
Idents(seu_obj_sub) <- seu_obj_sub$cluster_label_sub

marker_sheet_extended <- marker_tbl_extended %>%
  mutate(cluster_label_sub = cluster_label_extended[cluster_extended]) %>% 
  group_by(cluster_label_sub) %>%
  mutate(rank = row_number()) %>% 
  slice(1:50) %>% 
  select(rank, gene, cluster_label_sub) %>% 
  arrange(cluster_label_sub) %>% 
  spread(cluster_label_sub, gene)

cluster_n_tbl_full <- seu_obj$cluster_label_sub %>% 
  table() %>% 
  enframe("cluster_label_sub", "cluster_n") %>% 
  mutate(cluster_nrel = cluster_n/sum(cluster_n)) %>% 
  filter(cluster_label_sub != "NA")

marker_tbl_extended_annotated <- marker_tbl_extended %>% 
  mutate(cluster_label_sub = cluster_label_extended[cluster_extended]) %>%
  select(-cluster_extended) %>% 
  left_join(cluster_n_tbl_full, by = "cluster_label_sub") %>% 
  arrange(cluster_label_sub, desc(avg_logFC))

if (cell_sort == "CD45+") {
cells_to_keep <- colnames(seu_obj_sub)[!(str_detect(seu_obj_sub$cluster_label_sub, "dissociated|doublet|specific")) & !(str_detect(seu_obj_sub$cell_id, "CD45N"))]
}

if (cell_sort == "CD45-") {
cells_to_keep <- colnames(seu_obj_sub)[!(str_detect(seu_obj_sub$cluster_label_sub, "dissociated|doublet|specific")) & !(str_detect(seu_obj_sub$cell_id, "CD45P"))]
}


# seu_obj_sub_sub <- subset(seu_obj_sub, cells = cells_to_keep)
# seu_obj_sub_sub <- RunUMAP(seu_obj_sub_sub, dims = 1:50, reduction = "harmony", reduction.name = "umapharmony", reduction.key = "umapharmony_")
# write_rds(seu_obj_sub_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered_sub.rds"))

seu_obj_sub_sub <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered_sub.rds"))

seu_obj_sub_sub$cluster_label_sub <- setNames(cluster_label_extended[seu_obj_sub_sub$cluster_extended], seu_obj_sub_sub$cell_id)

write_rds(seu_obj_sub_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered_sub.rds"))

formattable::formattable(marker_sheet_extended)

rank	CD4.T.dysfunc.early	CD4.T.dysfunc.late.1	CD4.T.dysfunc.late.2	CD4.T.effector.memory	CD4.T.ISG	CD4.T.naive.centr.mem.1	CD4.T.naive.centr.mem.2	CD4.T.reg.1	CD4.T.reg.2	CD4.T.reg.3	CD4.T.reg.ISG	CD4.Th17.1	CD4.Th17.2	CD4.Th17.3	CD8.T.cytotoxic	CD8.T.dysfunc.early	CD8.T.dysfunc.ISG	CD8.T.dysfunc.late	CD8.T.effector.memory	CD8.T.ISG.early	CD8.T.ISG.late	CD8.T.naive.centr.mem	Cycling.CD4.T	Cycling.CD8.T.1	Cycling.CD8.T.2	Cycling.CD8.T.3	Cycling.CD8.T.4	Cycling.NK.1	Cycling.NK.2	Cycling.NK.3	dissociated_3	dissociated_4	doublet.Plasma.cell_1	gd.T.cell	NK.cytotoxic.GZMH	NK.cytotoxic.SPON2.1	NK.cytotoxic.SPON2.2	NK.reg.CCL3	NK.reg.CD56	NK.reg.CRTAM	NK.reg.IGFBP2	NK.reg.ISG	NK.reg.KRT81.KRT86.1	NK.reg.KRT81.KRT86.2	OV.009.specific
1	CXCL13	CXCL13	CXCL13	IL7R	IFIT3	CCR7	CCR7	RTKN2	TNFRSF4	TNFRSF4	ISG15	KLRB1	KLRB1	LST1	GZMH	CXCL13	CXCL13	CXCL13	GZMK	ISG15	IFIT2	CD8B	STMN1	STMN1	CXCL13	STMN1	STMN1	FCGR3A	CENPF	CENPF	MALAT1	EGR1	IGHG3	TRGC2	FGFBP2	FGFBP2	FGFBP2	CCL3	GNLY	XCL2	TYROBP	ISG15	GNLY	KRT81	IFITM3
2	IL6ST	PDCD1	NMB	CD40LG	ISG15	SELL	SELL	FOXP3	IL2RA	TNFRSF18	SAT1	IL7R	IL4I1	IL4I1	CD8A	HAVCR2	ISG15	KRT86	CD8B	IFIT3	IFIT3	IL7R	MKI67	HIST1H4C	HLA-DRA	TYMS	MKI67	IGFBP7	ASPM	ASPM	PTPRC	NR4A1	IGHG1	KLRC2	GZMH	SPON2	FCGR3A	AREG	AREG	GZMB	KLRC1	GNLY	FCER1G	KLRC1	IFIT1
3	FKBP5	CTLA4	TNFRSF4	ANXA1	MX2	EEF1B2	RPS8	UGP2	FOXP3	TNFRSF9	HERC5	CEBPD	IL7R	SCN1B	CD8B	MIR155HG	IFI6	GZMB	CD8A	IFIT1	IFIT1	CD8A	TUBA1B	TUBA1B	STMN1	PCLAF	TOP2A	KLRF1	MKI67	MKI67	TRMT112	KLRF1	IGHGP	XCL1	NKG7	FCGR3A	KLRF1	XCL1	KLRC1	CRTAM	FCER1G	IFIT3	KLRC1	TNFRSF18	IFIT3
4	NR3C1	MAF	CPM	GPR183	IFIT2	RPS12	RPS3A	TBC1D4	BATF	IL2RA	MX1	MYBL1	LTB	PCDH9	GZMK	PTMS	MX1	TNFRSF9	KLRG1	MX1	RSAD2	ANXA1	TUBB	MKI67	HLA-DRB1	DUT	ASPM	CX3CR1	TOP2A	TOP2A	C4orf3	FCGR3A	JCHAIN	ZNF683	CX3CR1	KLRF1	GNLY	XCL2	ITGAX	XCL1	TRDC	FCER1G	AREG	GNLY	RSAD2
5	MAF	ZBED2	IGFL2	KLF2	IFI44L	RPS6	RPS6	LTB	LTB	CTLA4	TNFRSF4	SLC4A10	TNFSF13B	SPINK2	CCL5	CD8A	LAG3	RBPJ	ITM2C	IFI6	MX1	TOB1	TOP2A	TUBB	FABP5	PCNA	NUSAP1	AKR1C3	KNL1	HMGB2	SF3B5	EGR2	IGLC2	LINC02446	KLF2	PRF1	PRF1	IER2	FCER1G	FABP5	AREG	MX1	TYROBP	FCER1G	ITGAD
6	CD40LG	IFNG	NR3C1	FOS	RSAD2	RPS3A	EEF1B2	IL2RA	TNFRSF18	BATF	CTLA4	CCR6	LST1	KIT	HLA-DRB1	LAG3	IFIT1	SRGAP3	PPP1R14B	MX2	ISG15	VIM	HMGB2	TOP2A	GAPDH	TUBA1B	TYMS	SPON2	NUSAP1	UBE2C	NSA2	FGFBP2	IGKC	HOPX	TRGC2	PTGDS	SPON2	FOS	XCL1	TNFRSF9	GNLY	IFIT1	TRDC	AREG	XCL1
7	AC004585.1	CD40LG	ICA1	LTB	IFIT1	RPS8	RPL34	IKZF2	CTLA4	GADD45A	SPATS2L	B3GALT2	CCR6	LTC4S	GZMA	GZMB	RSAD2	IFNG	DUSP2	ISG20	TNFSF10	YBX3	IL2RA	TYMS	IFNG	CLSPN	PCLAF	FGFBP2	CENPE	STMN1	SUMO1	SPON2	IGLC3	XCL2	PLEK	CX3CR1	KLRD1	EGR1	TYROBP	HSP90AB1	IGFBP2	TYROBP	KRT81	KRT86	FAM30A
8	TSHZ2	RBPJ	FKBP5	TPT1	MX1	LEF1	TCF7	PMAIP1	RTKN2	TNFRSF1B	IFI6	TPT1	SLC4A10	IL1R1	HLA-DPB1	CCL4L2	IFIT3	FAM3C	CST7	RSAD2	HERC5	S100A6	CENPF	NUSAP1	NASP	HELLS	UBE2C	GINS2	DLGAP5	TUBA1B	SEM1	IER3	IGHG2	CTSW	KLRG1	GNLY	TYROBP	CCL4	TRDC	PKM	KRT81	KLRC1	IGFBP2	RGS16	USP18
9	LIMS1	DUSP4	ITM2A	S100A4	TNFSF10	TCF7	RPS12	CTLA4	GADD45A	ICOS	IL2RA	RORA	CCL20	GSN	CXCR6	CXCR6	IFI44L	TIGIT	CMC1	IFIT2	MX2	FOS	NUSAP1	HMGB2	MCM5	MCM7	CLSPN	MLC1	GTSE1	UBE2S	SF3B6	CX3CR1	IGHG4	TRDC	ZEB2	GZMB	NKG7	NFKBIA	XCL2	NME1	KLRD1	IFITM3	TNFRSF18	XCL1	CMPK2
10	CD4	CD4	TSHZ2	NFKBIA	IL7R	RPS5	MAL	ARID5B	SOX4	FOXP3	FOXP3	ZBTB16	NCR3	LINC00299	HLA-DPA1	IFNG	SAMD9L	SPRY1	TRGC2	GZMK	ISG20	S100A10	PCLAF	HIST1H1B	COTL1	TUBB	CDK1	MYOM2	KIF14	TPX2	PSMD8	NFKBID	DERL3	KLRD1	PRF1	PLAC8	CLIC3	NR4A1	KLRD1	GAPDH	CLIC3	AREG	GSTP1	TYROBP	LINC02446
11	ITM2A	MIR155HG	IL6ST	TIMP1	IFI6	RPL32	RPL32	TIGIT	TNFRSF1B	PIM3	OAS1	DPP4	CTSH	AFF3	HLA-DRA	HLA-DRB1	PTMS	CCL3	CCL4	IFI44L	SAMD9L	ZFP36L2	HMGN2	CENPF	TUBA1B	NASP	BRCA1	E2F1	POLQ	CCNB1	GUK1	IGFBP7	MZB1	KIR3DL2	KLRD1	IGFBP7	GZMB	BTG2	CTSW	MIR155HG	CEBPD	RSAD2	XCL2	TRDC	ZNF683
12	CH25H	PRDM1	AHI1	JUN	STAT1	RPS18	KLF2	IL32	TIGIT	TYMP	IFIT1	KLRG1	AQP3	CD300LF	PTMS	RBPJ	OAS1	MIR155HG	CRTAM	MT2A	IFI6	MT-ND1	TYMS	ASPM	HELLS	MCM3	GTSE1	MCM7	TPX2	PTTG1	HNRNPF	CCL3	CD8A	TRGC1	PRSS23	FCER1G	PLAC8	DUSP2	CMC1	SRM	KLRB1	CD38	CTSW	ZNF683	IFIT2
13	TOX2	NR3C1	LIMS1	S100A10	HERC5	RPS13	RPL9	SAT1	LINC01943	LTB	IFIT3	JAML	DPP4	PLAT	ITM2C	CCL5	ISG20	PHLDA1	DTHD1	SAMD9L	MT2A	JUNB	H2AFZ	HMGN2	MCM4	FABP5	BIRC5	TYMS	BUB1B	NUSAP1	EIF1	IRF8	GZMK	CD7	FCGR3A	NKG7	CX3CR1	FOSB	NCAM1	RGS16	KRT86	IFI6	XCL1	HOPX	TRGC2
14	CYSLTR1	ITM2A	BIRC3	S100A11	ISG20	RPL13	RPS13	BATF	TBC1D4	TBC1D4	LTB	RPLP1	TPT1	DLL1	ZNF683	HLA-DRA	GZMB	CCND2	EOMES	OAS1	GBP1	MT-ND2	ASPM	PCLAF	RBPJ	MCM5	DLGAP5	S1PR5	HJURP	CENPE	ARPC1B	IER2	CD8B	NCR3	MYBL1	KLRD1	FCER1G	CD69	HOPX	NPW	XCL1	IFI44L	KRT86	CD7	IFI44L
15	PTPN13	TNFRSF18	CD200	CTSL	OAS1	RPL34	RPL21	AC005224.3	SAT1	IL1R2	TYMP	TMSB10	B3GALT2	LIF	COTL1	GZMH	CD8A	JAML	GZMM	EIF2AK2	OAS1	GZMK	PTTG1	H2AFZ	MCM6	GINS2	ORC6	MCM4	KIF11	KPNA2	UQCR11	CCL4	IGHA1	CD63	FLNA	CLIC3	PLEK	NR4A2	LINC00996	KDM6B	TXK	IFIT2	CD63	TMIGD2	OAS1
16	ST8SIA1	RGS1	GNG4	RPL34	OAS3	RPL21	RPL22	PBXIP1	UGP2	SAT1	RTKN2	RPS12	RORC	IL23R	CCL4	CCL3	LY6E	LAG3	LYAR	XAF1	NT5C3A	LMNA	CEP55	SMC4	PCNA	TK1	RRM2	RAMP1	SGO2	TUBB4B	UQCR10	AREG	CCL5	GNLY	EFHD2	TYROBP	EFHD2	NFKBID	IFITM3	TSPAN17	IL2RB	TRDC	MATK	KLRC2	HES4
17	IL6R	LIMS1	AC004585.1	JUNB	IFI44	NOSIP	RPL39	SELL	IKZF2	BIRC3	STAT1	TLE1	IL23R	TNFSF11	IFNG	PDCD1	EPSTI1	GAPDH	YBX3	STAT1	SAMD9	S100A4	TNFRSF4	RRM2	MIR155HG	MCM4	ZWINT	CHEK1	C21orf58	DLGAP5	SNU13	AKR1C3	ZNF66	IKZF2	ADGRG1	CCL3	CST7	KLRC1	KLRC2	GFOD1	MATK	TNFSF10	TMIGD2	CTSW	ISG15
18	SMCO4	LINC01943	CHN1	VIM	EIF2AK2	RPLP0	IL7R	BTG3	PMAIP1	PKM	ICOS	ERN1	LTK	RORC	CCL4L2	PHLDA1	HERC5	HAVCR2	CCL4L2	HERC5	IFI44L	CXCR4	DUT	HIST1H1D	MCM3	MKI67	RAD51AP1	CLSPN	KIF23	CKS2	SNRPD2	GNLY	PABPC1	TRG-AS1	CST7	AKR1C3	CCL3	TYROBP	IL2RB	PAICS	XCL2	CLIC3	GZMB	XCL2	OAS3
19	BTLA	ADAM19	TNFRSF18	RPLP0	GBP1	RPL5	RPL7	MIR4435-2HG	IL32	REL	EPSTI1	IL4I1	CA2	TGM2	HLA-DQA1	HLA-DPA1	HAVCR2	CTLA4	SH2D1A	LY6E	EPSTI1	TSC22D3	PCNA	UBE2C	HLA-DPA1	H2AFZ	CDT1	FCER1G	AC007240.1	HMMR	EIF3H	S1PR5	CCL4	CD52	GNLY	PLEK	ADGRG1	FCER1G	GZMB	MRTO4	NKG7	ISG20	CD7	SRGAP3	KLRC2
20	CPM	CD2	RNF19A	EEF1A1	SAT1	RPL9	RPL5	STAM	ICOS	IL1R1	TIGIT	EEF1A1	CEBPD	ZG16B	APOBEC3G	CD8B	OASL	CD8A	TC2N	IFI35	DDX58	LYAR	FOXP3	ATAD2	MKI67	ATAD2	KNL1	CHST2	KIF15	CCNB2	GTF3A	FCER1G	LTB	FCRL6	HOPX	EFHD2	ZEB2	IER3	TXK	EGR2	CD63	SAMD9L	IL2RB	CD63	MX2
21	PDE7B	TNFRSF4	MAF	RPS18	EPSTI1	RPS23	RPL37	BIRC3	TNFRSF9	ENTPD1	MX2	COTL1	RORA	IRAK3	CD3D	JAML	CD38	PTMS	CCL5	CD8A	OAS3	TNFRSF18	TPX2	TPX2	LAG3	HMGB2	HIST1H1B	STMN1	BRCA2	HMGN2	SEC61B	PRF1	DUSP23	HCST	C12orf75	ADGRG1	HOPX	TRDC	KLRF1	GEM	TMIGD2	BST2	CLIC3	CSF1	PLSCR1
22	GZMA	CCDC50	RBPJ	RPS12	MT2A	RPL3	RPS18	HPGD	LAYN	LINC01943	UGP2	HLA-DRB1	TNFRSF25	ENPP1	HLA-DQB1	HLA-DPB1	IFI35	ENTPD1	CXCR4	PLSCR1	PLSCR1	FOSB	ENO1	PCNA	TPI1	RANBP1	KIFC1	PLAC8	HMMR	SMC4	UXT	GZMB	CD3E	CCL5	C1orf21	CST7	AREG	KLRD1	CD7	CD72	CCL3	CEBPD	LAT2	CXXC5	MX1
23	CTSW	ETV7	SESN3	RPL39	SAMD9L	RPL7	RPL13	TNFRSF4	LINC02099	LTA	IFI44L	AQP3	RPLP1	PTGER3	CD74	VCAM1	PLSCR1	LINC01871	PPP2R5C	SAMD9	OASL	CD55	HIST1H1B	H2AFX	CLSPN	DEK	KIF11	GZMB	CIT	TUBB	NDUFA4	BTG2	GZMA	ZFP36L2	LITAF	CD247	TRDC	MAFF	SRGAP3	METTL1	CD7	GZMB	HOPX	SPRY1	IFI44
24	CCL5	MYO7A	BHLHE40-AS1	AQP3	CCR7	EEF1A1	RPL3	NABP1	CTSC	NAMPT	TBC1D4	TIGIT	TMIGD2	CA2	SH3BGRL3	LINC01871	IFI44	GOLIM4	GZMA	OASL	USP18	ZFP36	BATF	HMGB1	RANBP1	HNRNPAB	KIF15	TTC38	NCAPG	BIRC5	UQCRQ	KLRD1	CD3G	MATK	GZMB	AREG	CD247	KDM6B	CD300A	TNFRSF18	CTSW	OAS1	SRGAP3	CAPG	WARS
25	CCL4	TOX2	BATF	CD4	DDX58	RPS4X	RPS5	GADD45A	CORO1B	TIGIT	SAMD9	RPL13	TLE1	COL4A4	CD52	TNIP3	USP18	MYO1E	A2M-AS1	TRIM22	STAT1	DUSP1	CTLA4	HIST1H1C	ACTB	SMC2	HMMR	CTBP2	NCAPG2	CKS1B	RPL38	FOS	HLA-DPB1	IL32	PLAC8	HOPX	PTGDS	IGFBP2	MATK	GOLIM4	CMC1	PLSCR1	CXXC5	KIR2DL4	GMPR
26	NKG7	ITGA2	BTLA	TOB1	XAF1	CD8B	TSHZ2	CORO1B	ARID5B	PMAIP1	TNFRSF1B	NCR3	CD40LG	B3GALT5	CD3G	TNFRSF9	CXCR6	MYO7A	PDCD4	IRF7	IFI35	MT-CYB	SMC2	CKS1B	SNRPB	MCM6	PCNA	PLEK	CKAP2L	CDKN3	ATP5PD	NR4A3	PSMB1	LINC01871	ADRB2	MYOM2	KLRB1	CD160	SELL	NDFIP2	LAT2	KLRD1	KLRD1	LDLRAD4	IRF7
27	IGFBP4	CD200	CCDC50	TNFRSF25	TMEM123	RPL22	RPL11	IKZF4	CD27	F5	IKZF2	SPOCK2	ZBTB16	KIAA1211L	HLA-DRB5	HLA-DQA1	BST2	LAYN	PLSCR1	HELZ2	CMPK2	CD69	RRM2	DUT	ANXA5	SLBP	NCAPG	ADGRG1	ANLN	LGALS1	NDUFV2	GRASP	CD247	NUCB2	FCRL6	CHST2	C1orf21	EGR2	NKG7	PARK7	GSTP1	EPSTI1	ZNF683	LAT2	PADI4
28	CD200	BHLHE40-AS1	ZBED2	RPL11	CMPK2	RPS2	EEF1A1	AC133644.2	DUSP4	SOX4	TNFRSF18	LINC01871	JAML	HOXA5	THEMIS	CD63	STAT1	TNS3	IL6ST	OAS3	CD38	PPP1R14B	SMC4	TMPO	DEK	HMGN2	ATAD2	ORC6	KIF20B	ARL6IP1	UBE2L3	TYROBP	APOC1	IFITM2	ARL4C	PRSS23	CEBPD	IRF8	CD63	EIF5A	HOPX	IFI35	KLRC2	MATK	IFI6
29	CORO1B	PTPN13	CTLA4	FXYD5	SAMD9	RPL10A	RPL23	PHTF2	S100A4	SDC4	RSAD2	RPL32	S100A4	EMID1	LINC02446	TNFSF4	MT2A	CD63	IFI35	LAG3	EIF2AK2	AC020916.1	GAPDH	CLSPN	PTTG1	BRCA1	CENPU	MCM2	NUF2	CDC20	EIF3G	KDM6B	CRTAM	CRTAM	S1PR5	ZEB2	TXK	EGR3	KRT81	NOP16	SRGAP3	IRF7	NKG7	IL2RB	TRDC
30	CST7	PDE7B	CD4	RPS3A	DDX60L	RPL37	RPS23	ICOS	ENTPD1	DNPH1	SAMD9L	PDCD4	ERN1	NRP1	C12orf75	ID2	XAF1	SNX9	PKM	CMPK2	IFIH1	JUN	TK1	UBE2S	CD74	DNMT1	AURKB	PCNA	KIF4A	H2AFZ	ATP5F1C	PLAC8	HNRNPA1	ITGAX	SPON2	PTPN12	S1PR5	NR4A3	TNFRSF18	RGCC	CCL5	XCL2	SLC16A3	LINC01871	XCL2
31	FYB1	NKG7	TOX2	RPL3	TRIM22	RPL12	TPT1	DUSP16	BIRC3	CCND2	PMAIP1	HLA-DPB1	MYBL1	HOXA7	LAG3	LYST	TYMP	DGKH	EIF2AK2	GBP1	PARP14	TNFAIP3	PKM	KNL1	HLA-DQA1	CENPU	HELLS	CDT1	RPL10	KIF20B	ATP5F1B	CD160	CXCR6	KLRC4	KLF3	ABHD17A	IGFBP7	SPRY2	IGFBP2	DCTPP1	ITGA1	CCL3	NASP	ITGA1	IFI35
32	RUNX3	TSHZ2	METTL8	RPS6	HELZ2	RPL39	RPL30	CD27	CD4	LINC02099	USP18	ITGB1	SPOCK2	PDZK1	JAML	CD27	CCL5	AKAP5	TMEM123	IFI44	SAT1	MYBL1	PMAIP1	TUBB4B	TYMS	FEN1	MXD3	GNLY	ARHGAP11B	GTSE1	NDUFB4	NFKBIA	WNK1	CD9	S100A4	S1PR5	PRSS23	CD3D	CEBPD	ENO1	MCTP2	IL2RB	RGS16	PRMT9	ISG20
33	RNF19A	RNF19A	TNFRSF25	TNFAIP3	NT5C3A	LDHB	RPS28	PHACTR2	GLRX	TFRC	CD4	RPL34	FKBP11	NECTIN2	TRGC2	FABP5	GZMH	RGS2	ACP5	PARP14	IFI44	FCER1G	TMPO	SMC2	PTMS	CENPX	TPX2	PRF1	RPL13	TUBA1C	RARRES3	GADD45B	RGS1	PTGDR	PATL2	CTSW	PTPN12	CEBPD	IFITM2	IFNG	SH2D1B	HERC5	CD9	CTSA	HERC5
34	KLRD1	GAPDH	SMCO4	RPL9	CD40LG	RPSA	RPS20	CD4	MIR4435-2HG	DUSP4	SPOCK2	HLA-DPA1	EEF1A1	DACH1	TMSB4X	APOBEC3C	SAMD9	ITGAE	NDFIP2	EPSTI1	LAG3	GPR183	CLSPN	HIST1H1E	GINS2	DNAJC9	CDCA8	FEN1	RPS27	SGO2	SLC25A3	AC007952.4	ITM2C	MAF	TGFBR3	CCL4	GZMH	KRT81	GSTP1	TPI1	KLRF1	NT5C3A	STMN1	TXK	NT5C3A
35	TNFSF8	CORO1B	RILPL2	RPL32	USP18	RPL19	RPL38	RGS1	MAGEH1	ZBTB32	ARID5B	RPL3	KIF5C	SCARF1	LGALS1	TIGIT	EIF2AK2	DUSP4	NCR3	USP18	HELZ2	NCF1	TNFRSF18	BIRC5	DNMT1	GAPDH	TUBA1B	MCM5	RPL41	CDK1	EEF1D	MAFF	CD74	LYAR	CTSW	KLRB1	CMC1	ZBTB16	CAPN12	ZNF593	PRF1	PRF1	TXK	SLC16A3	CD7
36	SLC9A9	SAMSN1	PHACTR2	EEF1B2	LY6E	RPS16	RPLP0	FCMR	MAF	SH2D2A	NT5C3A	CD74	RUNX2	LINGO4	FABP5	ITGAE	GBP1	CLNK	COX5A	CD8B	PNPT1	KLF2	ATAD2	EZH2	MCM7	FAM111B	CENPW	HES6	RPS18	JPT1	TRIR	NR4A2	KIF13B	CD4	TTC38	TTC38	CTSW	GRASP	ZNF683	MALAT1	MAPK1	LY6E	ITGA1	NCR3	OAS2
37	PRF1	RPL3	THADA	RPL36	PARP14	S100B	LEF1	SPOCK2	CARD16	ID3	IFI44	IFNGR1	LINC01871	AL355974.2	CST7	HLA-DMA	HLA-DPA1	RGS1	ADGRG1	OAS2	XAF1	TYROBP	BIRC5	H2AFV	ACTG1	TPI1	CHEK1	ASCL2	STIL	KNL1	BLOC1S1	CCL4L2	C1QB	IFITM3	A2M-AS1	ITGB2	SAMD3	CCL3L1	KIR2DL4	PIM3	KLRC2	SAMD9	MCM5	NDFIP2	DDX58
38	PDCD1	CD3D	ARMH1	RPS8	IFIH1	RPS14	RPL36	CCNG2	SPOCK2	ARID5B	RNF213	RPS13	RPLP0	LPAR1	LGALS3	CRTAM	NT5C3A	NDFIP2	TMEM173	PARP9	TRIM22	CCL3	JPT1	CDK1	HAVCR2	CXCL13	SMC2	PTPN12	KIFC1	NUF2	COX5B	CD69	ARAP2	DGKD	TTC16	MYBL1	PTGDR	CD3G	CXXC5	ZBED2	CXXC5	XAF1	ADGRG3	CCL5	TNFSF10
39	CD8A	RPS12	CH25H	RPSA	IFI35	RPL11	RPS27	TRAC	IL1R1	CD4	SELL	RPLP0	TIGIT	PRAM1	SIT1	CTLA4	HLA-DRB1	TNIP3	CD247	DDX60	IRF7	CD4	UBE2C	KIF11	SMC2	SIVA1	FAM111B	ADAMTS1	RPL30	CEP55	SEC61G	PLEK	TPT1	CTLA4	ASCL2	TXK	TTC38	MAP3K8	ZFHX3	HSPA5	LINC00996	MX2	STARD3NL	CLNK	EPSTI1
40	CTLA4	RPL23A	AGFG1	LMNA	PLSCR1	LINC02446	PABPC1	ICA1	TYMP	CD27	PBXIP1	S100A4	HLA-DPB1	CD33	ANXA5	SNAP47	RBPJ	LINC01480	CSF1	RNF213	PPM1K	AHNAK	ANP32B	DEK	PFN1	CDT1	CENPF	CD300A	ECT2	KIF14	RPL31	NKG7	TNFRSF18	GPR183	FGR	XBP1	MYBL1	GADD45B	ID2	CD82	TNFRSF18	XCL1	KIR2DL4	CKLF	OASL
41	METTL8	ICOS	CYSLTR1	RPL13	IRF7	RPL23	RPS16	F5	PHACTR2	IKZF2	CD27	FYB1	RPS12	SLCO2A1	S100A6	RGS1	OAS3	NBL1	ALOX5AP	C19orf66	GBP4	MYADM	CENPE	DLGAP5	SNRPD1	SNRPB	KIF23	CDC45	RPS14	CENPA	OSTF1	CLIC3	SNTB1	KLRB1	LILRB1	PTGDR	AKR1C3	KLRF1	MCTP2	NOLC1	SAMD3	OASL	NCAM1	BCL2L11	ITGAX
42	ICA1	KLRD1	IGFBP4	RPL5	DDX60	RPL18A	RPS4X	FANK1	DNPH1	LAYN	OAS3	TNFSF13B	S100A6	AL034397.3	CD2	NKG7	IRF7	ID2	TPI1	PPM1K	DDX60L	KLF6	ICOS	ZWINT	SH3BGRL3	CENPM	UBE2T	CENPU	RPS12	MAD2L1	FTH1	RAP1GAP2	HLA-DPA1	CAPN12	VCL	CD300A	CEP78	SH2D1B	SPTSSB	CD160	CD247	USP18	MCTP2	CD247	DDX60L
43	TP53INP1	SRGN	SNX9	RPS27A	PPM1K	RPS28	RPLP2	S100A4	STAM	GBP2	NEAT1	HLA-DRA	IFNGR1	STAC	ITGA1	ITM2A	LINC02446	ZEB2	USP18	DDX58	DDX60	HAVCR2	CKS2	CENPM	CDCA7	TMEM106C	SGO1	GMNN	TPT1	CDCA8	ERH	MAP3K8	RPS13	CTSA	COL6A2	C1orf162	APMAP	ZFP36	MAFF	NPM3	STARD3NL	ERICH3	PRF1	CTSD	LGALS3BP
44	GZMH	IL6ST	ELMO1	RPL30	PNPT1	RPL4	SNHG8	NCF4	PHLDA1	ENO1	LINC02099	PLCB1	HLA-DRB1	TNFSF13B	GZMM	DUSP4	CD8B	RAB27A	CD226	SP110	RNF213	RGCC	TNFRSF1B	GTSE1	DNAJC9	E2F1	CDKN3	CCDC28B	RPLP1	HMGB3	NDUFS6	HOPX	IGHM	PRMT9	PLEKHG3	GZMH	CCL4	TXK	PTGDR	BZW2	PTGDR	EIF2AK2	SH2D1B	ADGRG3	BST2
45	TNFRSF4	RPL13	CD40LG	RPL38	OAS2	RPS27A	RPL12	AC084018.2	IL1R2	GCNT1	LINC01943	TRGC2	RPL34	TLE1	SLF1	ENTPD1	CMPK2	FABP5	ICOS	DDX60L	NUB1	PTGER4	KNL1	HIST2H2AC	CDT1	ATAD5	TUBB	TYROBP	RPL19	AURKA	SAP18	GSAP	GNLY	S100A6	CES1	LAIR2	MATK	CCL4L2	IL18	RELB	CD9	CD7	LINC00996	CD9	IFITM2
46	HCST	TYMP	MS4A6A	RPL37	SP110	RPLP2	RPS29	ZC2HC1A	SLAMF1	SYNGR2	LAYN	GZMB	PLCB1	FES	CRTAM	LINC02446	HELZ2	TNFRSF1B	ID2	BST2	XRN1	CTLA4	RANBP1	CDKN3	PDCD1	ZWINT	FEN1	KIR2DL3	RPL28	H2AFV	UQCRH	PPP1R15A	ZRANB2	CAPG	LINC02384	SH2D1B	ABHD17A	MATK	SH2D1B	MRPL12	CAPN12	OAS3	HELLS	GSTP1	LAMP3
47	GZMB	CTSW	LHFPL6	RPLP1	PARP9	RPL30	RPS27A	CCR4	F5	SERPINB9	ETV7	FKBP11	TMEM273	CSF2	CCR5	CLEC2D	RNF213	SNAP47	ANKRD28	IFI16	SP110	TNFRSF4	HELLS	FAM111B	CHEK1	NUDT1	NCAPG2	FGR	EEF1A1	PLK1	NOP10	CD3D	IKZF3	SCML4	COTL1	C1orf21	SH2D1B	LAT2	PLAC8	CCT2	NCAM1	WARS	ZBTB16	TRG-AS1	SPIRE1
48	CD84	SNX9	PDE7B	NOSIP	LGALS3BP	RPL29	RPS2	FTH1	LAIR2	PTP4A3	CMPK2	LTB	HLA-DPA1	BCAS1	APOBEC3C	AKAP5	VCAM1	ZBED2	STAT1	RARRES3	LAP3	DNAJB1	CKS1B	CENPE	NUDT1	DHFR	CENPE	LYN	RPS2	H2AFX	OAZ1	ZBTB16	UBE2D1	GLUL	CD27	ARL4C	MBP	NKG7	MMP25-AS1	RBPJ	APOBEC3G	IFI44	MCM4	KLRB1	SPINT2
49	ARID5B	SLAMF1	ST8SIA1	RPS16	LTB	MAL	RPL19	MALAT1	PHTF2	CCL20	F5	SATB1	KLRG1	RASSF8	HLA-DMA	CD3D	RPL13	ADAMTS6	TNFRSF25	SP100	IFI16	CD38	EZH2	ESCO2	SLBP	CHEK1	FANCI	PCLAF	RPS3	KIF11	RPL27A	DUSP2	TBCD	PRDM1	RGS1	KLF2	GK5	CD83	CAPG	RANBP1	METRNL	KRT81	TYMS	AFAP1L2	EIF2AK2
50	SESN3	RPS3	NPDC1	RPL21	PGAP1	TPT1	RPL35A	MT-CO2	AC005224.3	SOD1	IL1R1	S100A6	KLRD1	TOX2	GRAP2	CCND2	EEF1A1	BCL2L11	OAS3	CCL4	STAT2	TRDC	PHF19	SPC25	MYL6	RPA3	NUF2	EFHD2	RPS15A	SMC2	NDUFB7	CCL3L1	LARP7	GUK1	LTB	CD160	XBP1	REL	KLRB1	C1QBP	ZNF683	NKG7	MCM6	B3GNT7	SAMD9L

write_tsv(marker_sheet_extended, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_marker_sheet_full.tsv"))

write_tsv(marker_sheet_extended, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/supplementary_tables/", coi, "_marker_sheet_full.tsv"))

write_tsv(marker_tbl_extended_annotated, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_marker_table_annotated_full.tsv"))

write_tsv(marker_tbl_extended_annotated, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/supplementary_tables/", coi, "_marker_table_annotated_full.tsv"))

plot_data_sub_sub <- as_tibble(FetchData(seu_obj_sub_sub, c(myfeatures, "cluster_label_sub"))) %>% 
  left_join(meta_tbl, by = "sample") %>% 
  mutate(cell_id = colnames(seu_obj_sub_sub),
         cluster_label_sub = ordered(cluster_label_sub, levels = names(clrs$cluster_label_sub[[coi]])),
         ) %>%
  mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
         tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite)))
  
if (cell_sort == "CD45+") {
plot_data_sub_sub <- filter(plot_data_sub_sub, sort_parameters != "singlet, live, CD45-", !is.na(tumor_supersite))
}

if (cell_sort == "CD45-") {
plot_data_sub_sub <- filter(plot_data_sub_sub, sort_parameters != "singlet, live, CD45+", !is.na(tumor_supersite))
}

ggplot(plot_data_sub_sub, aes(umapharmony_1, umapharmony_2, color = cluster_label_sub)) + 
  common_layers_disc +
  # facet_wrap(~cluster_label_sub) +
  scale_color_manual(values = clrs$cluster_label_sub[[coi]]) +
  ggtitle("Sub cluster") 

4.1 add module scores and pathway scores

signature_modules_cd8 <- read_excel("_data/small/signatures/SPECTRUM v5 sub cluster markers.xlsx", sheet = 2, skip = 1, range = "M2:P100") %>% 
  gather(module, gene) %>% 
  na.omit() %>% 
  group_by(module) %>% 
  do(gene = c(.$gene)) %>% 
  {setNames(.$gene, .$module)} %>% 
  setNames(paste0(names(.), ".module"))

signature_modules_cin <- yaml::read_yaml("_data/small/signatures/cgas_sting.yaml") %>% 
  .[lapply(., length) > 1] %>% 
  setNames(paste0(names(.), ".module"))

# ## compute expression module scores
# for (i in 1:length(signature_modules_cd8)) {
#   seu_obj_sub_sub <- AddModuleScore(seu_obj_sub_sub, features = signature_modules_cd8[i], name = names(signature_modules_cd8)[i])
#   seu_obj_sub_sub[[names(signature_modules_cd8)[i]]] <- seu_obj_sub_sub[[paste0(names(signature_modules_cd8)[i], "1")]]
#   seu_obj_sub_sub[[paste0(names(signature_modules_cd8)[i], "1")]] <- NULL
#   print(paste(names(signature_modules_cd8)[i], "DONE"))
# }
# 
# for (i in 1:length(signature_modules_cin)) {
#   seu_obj_sub_sub <- AddModuleScore(seu_obj_sub_sub, features = signature_modules_cin[i], name = names(signature_modules_cin)[i])
#   seu_obj_sub_sub[[names(signature_modules_cin)[i]]] <- seu_obj_sub_sub[[paste0(names(signature_modules_cin)[i], "1")]]
#   seu_obj_sub_sub[[paste0(names(signature_modules_cin)[i], "1")]] <- NULL
#   print(paste(names(signature_modules_cin)[i], "DONE"))
# }
# 
# 
# ## compute progeny scores
# progeny_list <- seu_obj_sub_sub@assays$RNA@data[VariableFeatures(seu_obj_sub_sub),] %>%
#   as.matrix %>%
#   progeny %>%
#   as.data.frame %>%
#   as.list
# 
# names(progeny_list) <- make.names(paste0(names(progeny_list), ".pathway"))
# 
# for (i in 1:length(progeny_list)) {
#   seu_obj_sub_sub <- AddMetaData(seu_obj_sub_sub,
#                                  metadata = progeny_list[[i]],
#                                  col.name = names(progeny_list)[i])
# }
# 
# write_rds(seu_obj_sub_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered_annotated.rds"))

seu_obj_sub_sub <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered_annotated.rds"))

seu_obj_sub_sub$cluster_label_sub <- setNames(cluster_label_extended[seu_obj_sub_sub$cluster_extended], seu_obj_sub_sub$cell_id)

write_rds(seu_obj_sub_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_processed_filtered_annotated.rds"))

4.2 marker heatmap

marker_top_tbl <- marker_sheet_extended[,-1] %>% 
  slice(1:10) %>% 
  as.list %>% 
  .[!str_detect(names(.), "doublet|dissociated")] %>% 
  enframe("cluster_label_x", "gene") %>% 
  unnest(gene)

plot_data_markers <- as_tibble(FetchData(seu_obj_sub_sub, c("cluster_label_sub", myfeatures, unique(marker_top_tbl$gene)))) %>% 
  gather(gene, value, -c(1:(length(myfeatures)+1))) %>% 
  left_join(meta_tbl, by = "sample") %>% 
  mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
         tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite))) %>% 
  mutate(cluster_label_sub = ordered(cluster_label_sub, levels = names(clrs$cluster_label_sub[[coi]]))) %>% 
  group_by(cluster_label_sub, gene) %>% 
  summarise(value = mean(value, na.rm = T)) %>% 
  group_by(gene) %>% 
  mutate(value = scales::rescale(value)) %>% 
  left_join(marker_top_tbl, by = "gene") %>% 
  mutate(cluster_label_x = ordered(cluster_label_x, levels = rev(names(clrs$cluster_label_sub[[coi]]))))

ggplot(plot_data_markers) +
  geom_tile(aes(gene, cluster_label_sub, fill = value)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
  facet_wrap(~cluster_label_x, scales = "free") +
  scale_fill_gradientn(colors = viridis(9)) +
  labs(fill = "Scaled\nexpression") +
  theme(aspect.ratio = 1,
        axis.line = element_blank(),
        axis.ticks = element_blank(),
        axis.title = element_blank())

# ggsave(paste0("_fig/002_marker_heatmap_", coi, ".pdf"), width = nrow(marker_top_tbl)/6, height = 5)

4.3 composition

4.3.1 per site

comp_site_tbl <- plot_data_sub_sub %>%
  filter(!is.na(tumor_supersite)) %>% 
  group_by(cluster_label_sub, tumor_supersite) %>%
  tally %>%
  group_by(tumor_supersite) %>%
  mutate(nrel = n/sum(n)*100) %>%
  ungroup

pnrel_site <- ggplot(comp_site_tbl) +
  geom_bar(aes(tumor_supersite, nrel, fill = cluster_label_sub),
           stat = "identity", position = position_stack()) +
  scale_y_continuous(expand = c(0, 0)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
  labs(fill = "Cluster", y = "Fraction [%]", x = "") +
  scale_fill_manual(values = clrs$cluster_label_sub[[coi]])

pnabs_site <- ggplot(comp_site_tbl) +
  geom_bar(aes(tumor_supersite, n, fill = cluster_label_sub),
           stat = "identity", position = position_stack()) +
  scale_y_continuous(expand = c(0, 0)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
  labs(fill = "Cluster", y = "# cells", x = "") +
  scale_fill_manual(values = clrs$cluster_label_sub[[coi]])

plot_grid(pnabs_site, pnrel_site, ncol = 2, align = "h")

# ggsave(paste0("_fig/02_deep_dive_", coi, "_comp_site.pdf"), width = 8, height = 4)

4.3.2 per sample

comp_tbl_sample_sort <- plot_data_sub_sub %>% 
  group_by(tumor_subsite, tumor_supersite, patient_id, consensus_signature, therapy, cluster_label_sub) %>% 
  tally %>% 
  group_by(tumor_subsite, tumor_supersite, patient_id, consensus_signature, therapy) %>% 
  mutate(nrel = n/sum(n)*100,
         nsum = sum(n),
         log10n = log10(n),
         label_supersite = "Site",
         label_mutsig = "Signature",
         label_therapy = "Rx") %>% 
  ungroup %>% 
  arrange(desc(therapy), tumor_supersite) %>% 
  mutate(tumor_subsite_rx = paste0(tumor_subsite, "_", therapy)) %>% 
  mutate(tumor_subsite = ordered(tumor_subsite, levels = unique(tumor_subsite)),
         tumor_subsite_rx = ordered(tumor_subsite_rx, levels = unique(tumor_subsite_rx))) %>% 
  arrange(patient_id) %>% 
  mutate(label_patient_id = ifelse(as.logical(as.numeric(fct_inorder(as.character(patient_id)))%%2), "Patient1", "Patient2"))

sample_id_x_tbl <- plot_data_sub %>% 
  mutate(sort_short_x = cell_sort) %>% 
  distinct(patient_id, sort_short_x, tumor_subsite, therapy, sample) %>% 
  unite(sample_id_x, patient_id, sort_short_x, tumor_subsite, therapy) %>% 
  arrange(sample_id_x)

comp_tbl_sample_sort %>% 
  mutate(sort_short_x = cell_sort) %>% 
  unite(sample_id_x, patient_id, sort_short_x, tumor_subsite_rx) %>% 
  select(sample_id_x, cluster_label_sub, n, nrel, nsum) %>% 
  left_join(sample_id_x_tbl, by = "sample_id_x") %>% 
  write_tsv(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_subtype_compositions.tsv"))

ybreaks <- c(500, 1000, 2000, 4000, 6000, 8000, 10000, 15000, 20000)

max_cells_per_sample <- max(comp_tbl_sample_sort$nsum)
ymaxn <- ybreaks[max_cells_per_sample < ybreaks][1]

comp_plot_wrapper <- function(y = "nrel", switch = NULL) {
  if (y == "nrel") ylab <- paste0("Fraction\nof cells [%]")
  if (y == "n") ylab <- paste0("Number\nof cells")
  p <- ggplot(comp_tbl_sample_sort, 
              aes_string("tumor_subsite_rx", y, fill = "cluster_label_sub")) + 
    facet_grid(~patient_id, space = "free", scales = "free", switch = switch) +
    coord_cartesian(clip = "off") + 
    scale_fill_manual(values = clrs$cluster_label_sub[[coi]]) + 
    theme(axis.text.x = element_blank(),
          axis.title.y = element_text(angle = 0, vjust = 0.5, hjust = 0.5, 
                                      margin = margin(0, -0.4, 0, 0, unit = "npc")),
          axis.ticks.x = element_blank(),
          axis.title.x = element_blank(),
          axis.line.x = element_blank(),
          strip.text.y = element_blank(),
          strip.text.x = element_blank(),
          strip.background.y = element_blank(),
          strip.background.x = element_blank(),
          plot.margin = margin(0, 0, 0, 0, "npc")) + 
    labs(x = "",
         y = ylab) +
    guides(fill = FALSE)
  if (y == "nrel") p <- p + 
    geom_bar(stat = "identity") +
    scale_y_continuous(expand = c(0, 0), 
                       breaks = c(0, 50, 100), 
                       labels = c("0", "50", "100"))
  if (y == "n") p <- p + 
    geom_bar(stat = "identity", position = position_stack()) +
    scale_y_continuous(expand = c(0, 0), 
                       breaks = c(0, ymaxn/2, ymaxn),
                       limits = c(0, ymaxn),
                       labels = c("", ymaxn/2, ymaxn)) +
    expand_limits(y = c(0, ymaxn)) +
    theme(panel.grid.major.y = element_line(linetype = 1, color = "grey90", size = 0.5))
  return(p)
} 

common_label_layers <- list(
  geom_tile(color = "white", size = 0.15),
  theme(axis.text.x = element_blank(),
        axis.ticks = element_blank(),
        axis.title.x = element_blank(),
        axis.line.x = element_blank(),
        strip.text = element_blank(),
        strip.background = element_blank(),
        plot.margin = margin(0, 0, 0, 0, "npc")),
  scale_y_discrete(expand = c(0, 0)),
  labs(y = ""),
  guides(fill = FALSE),
  facet_grid(~patient_id, 
             space = "free", scales = "free")
)

comp_label_site <- ggplot(distinct(comp_tbl_sample_sort, tumor_subsite_rx, therapy, tumor_supersite, label_supersite, patient_id), 
       aes(tumor_subsite_rx, label_supersite, 
           fill = tumor_supersite)) + 
  scale_fill_manual(values = clrs$tumor_supersite) +
  common_label_layers

comp_label_rx <- ggplot(distinct(comp_tbl_sample_sort, tumor_subsite_rx, therapy, tumor_supersite, label_therapy, consensus_signature, patient_id), 
       aes(tumor_subsite_rx, label_therapy, 
           fill = therapy)) + 
  scale_fill_manual(values = c(`post-Rx` = "gold3", `pre-Rx` = "steelblue")) +
  common_label_layers

comp_label_mutsig <- ggplot(distinct(comp_tbl_sample_sort, tumor_subsite_rx, therapy, tumor_supersite, label_mutsig, consensus_signature, patient_id), 
       aes(tumor_subsite_rx, label_mutsig, 
           fill = consensus_signature)) + 
  scale_fill_manual(values = clrs$consensus_signature) +
  common_label_layers

patient_label_tbl <- distinct(comp_tbl_sample_sort, patient_id, .keep_all = T)

comp_label_patient_id <- ggplot(patient_label_tbl, aes(tumor_subsite_rx, label_patient_id)) + 
  scale_fill_manual(values = clrs$consensus_signature) +
  geom_text(aes(tumor_subsite_rx, label_patient_id, label = patient_id)) +
  facet_grid(~patient_id, 
             space = "free", scales = "free") +
  coord_cartesian(clip = "off") + 
  theme_void() +
  theme(strip.text = element_blank(),
        strip.background = element_blank(),
        plot.margin = margin(0, 0, 0, 0, "npc"))

hist_plot_wrapper <- function(x = "nrel") {
  if (x == "nrel") {
    xlab <- paste0("Fraction of cells [%]")
    bw <- 5
  }
  if (x == "log10n") {
    xlab <- paste0("Number of cells")
    bw <- 0.2
  }
  p <- ggplot(comp_tbl_sample_sort) +
    ggridges::geom_density_ridges(
      aes_string(x, "cluster_label_sub", fill = "cluster_label_sub"), color = "black",
      stat = "binline", binwidth = bw, scale = 3) +
    facet_grid(label_supersite~., 
               space = "free", scales = "free") +
    scale_fill_manual(values = clrs$cluster_label_sub[[coi]]) + 
    theme(axis.text.y = element_blank(),
          axis.ticks.y = element_blank(),
          axis.title.y = element_blank(),
          axis.line.y = element_blank(),
          strip.text = element_blank(),
          strip.background = element_blank(),
          plot.margin = margin(0, 0, 0, 0, "npc")) +
    scale_x_continuous(expand = c(0, 0)) +
    scale_y_discrete(expand = c(0, 0)) +
    guides(fill = FALSE) +
    labs(x = xlab)
  if (x == "log10n") p <- p + expand_limits(x = c(0, 3)) + 
    scale_x_continuous(expand = c(0, 0), 
                       labels = function(x) parse(text = paste("10^", x)))
  return(p)
}

pcomp1 <- comp_plot_wrapper("n")
pcomp2 <- comp_plot_wrapper("nrel")

pcomp_grid <- plot_grid(comp_label_patient_id, 
                        pcomp1, pcomp2, 
                        comp_label_site, comp_label_rx, comp_label_mutsig,
                        ncol = 1, align = "v", 
                        rel_heights = c(0.15, 0.33, 0.33, 0.06, 0.06, 0.06))

phist1 <- hist_plot_wrapper("log10n")

pcomp_hist_grid <- ggdraw() +
  draw_plot(pcomp_grid, x = 0.01, y = 0, width = 0.85, height = 1) +
  draw_plot(phist1, x = 0.87, y = 0.05, width = 0.12, height = 0.8)

pcomp_hist_grid

# ggsave(paste0("_fig/02_composition_v7_",coi,".pdf"), pcomp_hist_grid, width = 10, height = 2)

4.3.3 site specific cluster enrichment

comp_tbl_z <- comp_tbl_sample_sort %>% 
  filter(therapy == "pre-Rx",
         !(tumor_supersite %in% c("Ascites", "Other"))) %>% 
  group_by(patient_id, cluster_label_sub) %>% 
  arrange(patient_id, cluster_label_sub, nrel) %>% 
  mutate(rank = row_number(nrel),
         z_rank = scales::rescale(rank)) %>% 
  mutate(mean_nrel = mean(nrel, na.rm = T),
         sd_nrel = sd(nrel, na.rm = T),
         z_nrel = (nrel - mean_nrel) / sd_nrel) %>% 
  ungroup()

ggplot(comp_tbl_z) +
  geom_boxplot(aes(tumor_supersite, z_nrel, color = tumor_supersite), 
               outlier.shape = NA) +
  geom_point(aes(tumor_supersite, z_nrel, color = tumor_supersite), position = position_jitter(), size = 0.1) +
  facet_grid(~cluster_label_sub, scales = "free_x", space = "free_x") +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5),
        aspect.ratio = 5) +
  scale_color_manual(values = clrs$tumor_supersite)

ggplot(comp_tbl_z) +
  geom_boxplot(aes(tumor_supersite, z_rank, color = tumor_supersite), 
               outlier.shape = NA) +
  geom_point(aes(tumor_supersite, z_rank, color = tumor_supersite), position = position_jitter(), size = 0.1) +
  facet_grid(~cluster_label_sub, scales = "free_x", space = "free_x") +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5),
        aspect.ratio = 5) +
  scale_color_manual(values = clrs$tumor_supersite)

4.4 CD8 T cell diffusion maps

# subset CD8 cells
seu_obj_cd8 <- seu_obj_sub_sub %>%
  subset(subset = cluster_label_sub %in% grep("^CD8.T", unique(seu_obj_sub_sub$cluster_label_sub), value = T)) %>%
  RunUMAP(dims = 1:50, reduction = "harmony", reduction.name = "umapharmony", reduction.key = "umapharmony_")

write_rds(seu_obj_cd8, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/CD8.T_processed.rds")

seu_obj_cd8 <- read_rds("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/CD8.T_processed.rds")

seu_obj_cd8$cluster_label_sub <- setNames(cluster_label_extended[seu_obj_cd8$cluster_extended], seu_obj_cd8$cell_id)

write_rds(seu_obj_cd8, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/CD8.T_processed.rds")

root_cell <- "SPECTRUM-OV-045_S1_CD45P_LEFT_OVARY_TCATTCAGTAGTCGGA"
seu_obj_cd8_sub <- subset(seu_obj_cd8, cells = c(root_cell, colnames(seu_obj_cd8)[colnames(seu_obj_cd8)!=root_cell][-1]))

dc_obj <- DiffusionMap(seu_obj_cd8_sub@reductions$harmony@cell.embeddings, k = 100)
dc_mat <- dc_obj@eigenvectors
colnames(dc_mat) <- paste0("DC_", 1:ncol(dc_mat))
seu_obj_cd8_sub[["DC"]] <- CreateDimReducObject(embeddings = dc_mat, key = "DC_", assay = DefaultAssay(seu_obj_cd8_sub))

write_rds(seu_obj_cd8_sub, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/CD8.T_processed_filtered.rds")

dpt_obj <- DPT(dc_obj)
write_rds(dpt_obj, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/CD8.T_sub_dpt.rds")
seu_obj_cd8_sub$DPT1 <- dpt_obj$DPT1
write_rds(seu_obj_cd8_sub, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/CD8.T_processed_filtered.rds")

seu_obj_cd8_sub <- read_rds("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/CD8.T_processed_filtered.rds")

seu_obj_cd8_sub$cluster_label_sub <- setNames(cluster_label_extended[seu_obj_cd8_sub$cluster_extended], seu_obj_cd8_sub$cell_id)

write_rds(seu_obj_cd8_sub, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/CD8.T_processed_filtered.rds")

# plot_data_sub_sub <- as_tibble(FetchData(seu_obj_sub_sub, c(myfeatures, "cluster_label"))) %>% 
#   left_join(meta_tbl, by = "sample") %>% 
#   mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
#          tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite))) %>% 
#   mutate(cell_id = colnames(seu_obj_sub_sub),
#          cluster_label = ordered(cluster_label, levels = my_subtypes),
#          )
#   
# if (cell_sort == "CD45+") {
# plot_data_sub_sub <- filter(plot_data_sub_sub, sort_parameters != "singlet, live, CD45-", !is.na(tumor_supersite))
# }
# 
# if (cell_sort == "CD45-") {
# plot_data_sub_sub <- filter(plot_data_sub_sub, sort_parameters != "singlet, live, CD45+", !is.na(tumor_supersite))
# }
# 
# ggplot(plot_data_sub_sub, aes(umapharmony_1, umapharmony_2, color = cluster_label)) + 
#   common_layers_disc +
#   ggtitle("Sub cluster") +
#   #facet_wrap(~cluster_label) +
#   scale_color_manual(values = clrs$cluster_label[[coi]])
# 
# ggplot(plot_data_sub_sub, aes(umapharmony_1, umapharmony_2, color = patient_id_short)) + 
#   common_layers_disc +
#   ggtitle("Patient") +
#   # facet_wrap(~therapy) +
#   scale_color_manual(values = clrs$patient_id_short)
# 
# ggplot(plot_data_sub_sub, aes(umapharmony_1, umapharmony_2, color = tumor_supersite)) + 
#   # geom_point(aes(umapharmony_1, umapharmony_2), 
#   #            color = "grey90", size = 0.01, 
#   #            data = select(plot_data_sub_sub, -tumor_supersite)) +
#   common_layers_disc +
#   ggtitle("Site") +
#   # facet_wrap(~therapy) +
#   scale_color_manual(values = clrs$tumor_supersite)
# 
# write_tsv(select(plot_data_sub_sub, cell_id, everything(), -umapharmony_1, -umapharmony_2, -contains("RNA_")), paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v7/", coi, "_embedding_sub.tsv"))

4.4.1 UMAP & DC

plot_data_cd8_sub <- as_tibble(FetchData(seu_obj_cd8_sub, c(myfeatures, "cluster_label_sub", "DC_1", "DC_2"))) %>% 
  left_join(meta_tbl, by = "sample") %>% 
  mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
         tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite))) %>% 
  mutate(cell_id = colnames(seu_obj_cd8_sub),
         cluster_label_sub = ordered(cluster_label_sub, levels = names(clrs$cluster_label_sub[[coi]])))
  
if (cell_sort == "CD45+") {
plot_data_cd8_sub <- filter(plot_data_cd8_sub, sort_parameters != "singlet, live, CD45-", !is.na(tumor_supersite))
}

if (cell_sort == "CD45-") {
plot_data_cd8_sub <- filter(plot_data_cd8_sub, sort_parameters != "singlet, live, CD45+", !is.na(tumor_supersite))
}

ggplot(plot_data_cd8_sub, aes(umapharmony_1, umapharmony_2, color = cluster_label_sub)) + 
  common_layers_disc +
  ggtitle("Sub cluster") +
  #facet_wrap(~cluster_label_sub) +
  scale_color_manual(values = clrs$cluster_label_sub[[coi]])

ggplot(plot_data_cd8_sub, aes(umapharmony_1, umapharmony_2, color = patient_id_short)) + 
  common_layers_disc +
  ggtitle("Patient") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$patient_id_short)

ggplot(plot_data_cd8_sub, aes(umapharmony_1, umapharmony_2, color = tumor_supersite)) + 
  common_layers_disc +
  ggtitle("Site") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$tumor_supersite)

ggplot(plot_data_cd8_sub, aes(DC_1, DC_2, color = cluster_label_sub)) + 
  common_layers_disc +
  ggtitle("Sub cluster") +
  #facet_wrap(~cluster_label_sub) +
  scale_color_manual(values = clrs$cluster_label_sub[[coi]])

ggplot(plot_data_cd8_sub, aes(DC_1, DC_2, color = patient_id_short)) + 
  common_layers_disc +
  ggtitle("Patient") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$patient_id_short)

ggplot(plot_data_cd8_sub, aes(DC_1, DC_2, color = tumor_supersite)) + 
  common_layers_disc +
  ggtitle("Site") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$tumor_supersite)

5 session info

devtools::session_info()

## ─ Session info ───────────────────────────────────────────────────────────────
##  setting  value                       
##  version  R version 3.6.2 (2019-12-12)
##  os       Debian GNU/Linux 10 (buster)
##  system   x86_64, linux-gnu           
##  ui       X11                         
##  language (EN)                        
##  collate  en_US.UTF-8                 
##  ctype    en_US.UTF-8                 
##  tz       Etc/UTC                     
##  date     2021-03-08                  
## 
## ─ Packages ───────────────────────────────────────────────────────────────────
##  package              * version    date       lib
##  abind                  1.4-5      2016-07-21 [2]
##  ape                    5.3        2019-03-17 [2]
##  assertthat             0.2.1      2019-03-21 [2]
##  backports              1.1.10     2020-09-15 [1]
##  bibtex                 0.4.2.2    2020-01-02 [2]
##  Biobase                2.46.0     2019-10-29 [2]
##  BiocGenerics           0.32.0     2019-10-29 [2]
##  BiocParallel           1.20.1     2019-12-21 [2]
##  bitops                 1.0-6      2013-08-17 [2]
##  boot                   1.3-24     2019-12-20 [3]
##  broom                  0.7.2      2020-10-20 [1]
##  callr                  3.4.2      2020-02-12 [1]
##  car                    3.0-8      2020-05-21 [1]
##  carData                3.0-4      2020-05-22 [1]
##  caTools                1.17.1.4   2020-01-13 [2]
##  cellranger             1.1.0      2016-07-27 [2]
##  class                  7.3-15     2019-01-01 [3]
##  cli                    2.0.2      2020-02-28 [1]
##  cluster                2.1.0      2019-06-19 [3]
##  codetools              0.2-16     2018-12-24 [3]
##  colorblindr          * 0.1.0      2020-01-13 [2]
##  colorspace           * 1.4-2      2019-12-29 [2]
##  cowplot              * 1.0.0      2019-07-11 [2]
##  crayon                 1.3.4      2017-09-16 [1]
##  curl                   4.3        2019-12-02 [2]
##  data.table             1.12.8     2019-12-09 [2]
##  DBI                    1.1.0      2019-12-15 [2]
##  dbplyr                 2.0.0      2020-11-03 [1]
##  DelayedArray           0.12.2     2020-01-06 [2]
##  DEoptimR               1.0-8      2016-11-19 [1]
##  desc                   1.2.0      2018-05-01 [2]
##  destiny              * 3.0.1      2020-01-16 [1]
##  devtools               2.2.1      2019-09-24 [2]
##  digest                 0.6.25     2020-02-23 [1]
##  dplyr                * 1.0.2      2020-08-18 [1]
##  e1071                  1.7-3      2019-11-26 [1]
##  ellipsis               0.3.1      2020-05-15 [1]
##  evaluate               0.14       2019-05-28 [2]
##  fansi                  0.4.1      2020-01-08 [2]
##  farver                 2.0.3      2020-01-16 [1]
##  fitdistrplus           1.0-14     2019-01-23 [2]
##  forcats              * 0.5.0      2020-03-01 [1]
##  foreign                0.8-74     2019-12-26 [3]
##  formattable            0.2.0.1    2016-08-05 [1]
##  fs                     1.5.0      2020-07-31 [1]
##  future                 1.15.1     2019-11-25 [2]
##  future.apply           1.4.0      2020-01-07 [2]
##  gbRd                   0.4-11     2012-10-01 [2]
##  gdata                  2.18.0     2017-06-06 [2]
##  generics               0.0.2      2018-11-29 [2]
##  GenomeInfoDb           1.22.0     2019-10-29 [2]
##  GenomeInfoDbData       1.2.2      2020-01-14 [2]
##  GenomicRanges          1.38.0     2019-10-29 [2]
##  ggplot.multistats      1.0.0      2019-10-28 [1]
##  ggplot2              * 3.3.2      2020-06-19 [1]
##  ggrepel                0.8.1      2019-05-07 [2]
##  ggridges               0.5.2      2020-01-12 [2]
##  ggthemes               4.2.0      2019-05-13 [1]
##  globals                0.12.5     2019-12-07 [2]
##  glue                   1.3.2      2020-03-12 [1]
##  gplots                 3.0.1.2    2020-01-11 [2]
##  gridExtra              2.3        2017-09-09 [2]
##  gtable                 0.3.0      2019-03-25 [2]
##  gtools                 3.8.1      2018-06-26 [2]
##  haven                  2.3.1      2020-06-01 [1]
##  hexbin                 1.28.0     2019-11-11 [2]
##  hms                    0.5.3      2020-01-08 [1]
##  htmltools              0.5.1.1    2021-01-22 [1]
##  htmlwidgets            1.5.1      2019-10-08 [2]
##  httr                   1.4.2      2020-07-20 [1]
##  ica                    1.0-2      2018-05-24 [2]
##  igraph                 1.2.5      2020-03-19 [1]
##  IRanges                2.20.2     2020-01-13 [2]
##  irlba                  2.3.3      2019-02-05 [2]
##  jsonlite               1.7.1      2020-09-07 [1]
##  KernSmooth             2.23-16    2019-10-15 [3]
##  knitr                  1.26       2019-11-12 [2]
##  knn.covertree          1.0        2019-10-28 [1]
##  labeling               0.3        2014-08-23 [2]
##  laeken                 0.5.1      2020-02-05 [1]
##  lattice                0.20-38    2018-11-04 [3]
##  lazyeval               0.2.2      2019-03-15 [2]
##  leiden                 0.3.1      2019-07-23 [2]
##  lifecycle              0.2.0      2020-03-06 [1]
##  listenv                0.8.0      2019-12-05 [2]
##  lmtest                 0.9-37     2019-04-30 [2]
##  lsei                   1.2-0      2017-10-23 [2]
##  lubridate              1.7.9.2    2020-11-13 [1]
##  magrittr             * 2.0.1      2020-11-17 [1]
##  MASS                   7.3-51.5   2019-12-20 [3]
##  Matrix                 1.2-18     2019-11-27 [3]
##  matrixStats            0.56.0     2020-03-13 [1]
##  memoise                1.1.0      2017-04-21 [2]
##  metap                  1.2        2019-12-08 [2]
##  mnormt                 1.5-5      2016-10-15 [2]
##  modelr                 0.1.8      2020-05-19 [1]
##  multcomp               1.4-12     2020-01-10 [2]
##  multtest               2.42.0     2019-10-29 [2]
##  munsell                0.5.0      2018-06-12 [2]
##  mutoss                 0.1-12     2017-12-04 [2]
##  mvtnorm                1.0-12     2020-01-09 [2]
##  nlme                   3.1-143    2019-12-10 [3]
##  nnet                   7.3-12     2016-02-02 [3]
##  npsurv                 0.4-0      2017-10-14 [2]
##  numDeriv               2016.8-1.1 2019-06-06 [2]
##  openxlsx               4.1.5      2020-05-06 [1]
##  pbapply                1.4-2      2019-08-31 [2]
##  pcaMethods             1.78.0     2019-10-29 [2]
##  pillar                 1.4.6      2020-07-10 [1]
##  pkgbuild               1.0.6      2019-10-09 [2]
##  pkgconfig              2.0.3      2019-09-22 [1]
##  pkgload                1.0.2      2018-10-29 [2]
##  plotly                 4.9.1      2019-11-07 [2]
##  plotrix                3.7-7      2019-12-05 [2]
##  plyr                   1.8.5      2019-12-10 [2]
##  png                    0.1-7      2013-12-03 [2]
##  prettyunits            1.1.1      2020-01-24 [1]
##  processx               3.4.2      2020-02-09 [1]
##  progeny              * 1.11.3     2020-10-22 [1]
##  proxy                  0.4-24     2020-04-25 [1]
##  ps                     1.3.2      2020-02-13 [1]
##  purrr                * 0.3.4      2020-04-17 [1]
##  R.methodsS3            1.7.1      2016-02-16 [2]
##  R.oo                   1.23.0     2019-11-03 [2]
##  R.utils                2.9.2      2019-12-08 [2]
##  R6                     2.4.1      2019-11-12 [1]
##  ranger                 0.12.1     2020-01-10 [1]
##  RANN                   2.6.1      2019-01-08 [2]
##  rappdirs               0.3.1      2016-03-28 [2]
##  RColorBrewer           1.1-2      2014-12-07 [2]
##  Rcpp                   1.0.4      2020-03-17 [1]
##  RcppAnnoy              0.0.16     2020-03-08 [1]
##  RcppEigen              0.3.3.7.0  2019-11-16 [2]
##  RcppHNSW               0.2.0      2019-09-20 [2]
##  RcppParallel           4.4.4      2019-09-27 [2]
##  RCurl                  1.98-1.1   2020-01-19 [1]
##  Rdpack                 0.11-1     2019-12-14 [2]
##  readr                * 1.4.0      2020-10-05 [1]
##  readxl               * 1.3.1      2019-03-13 [2]
##  rematch                1.0.1      2016-04-21 [2]
##  remotes                2.1.0      2019-06-24 [2]
##  reprex                 0.3.0      2019-05-16 [2]
##  reshape2               1.4.3      2017-12-11 [2]
##  reticulate             1.14       2019-12-17 [2]
##  rio                    0.5.16     2018-11-26 [1]
##  rlang                  0.4.8      2020-10-08 [1]
##  rmarkdown              2.0        2019-12-12 [2]
##  robustbase             0.93-6     2020-03-23 [1]
##  ROCR                   1.0-7      2015-03-26 [2]
##  rprojroot              1.3-2      2018-01-03 [2]
##  RSpectra               0.16-0     2019-12-01 [2]
##  rstudioapi             0.11       2020-02-07 [1]
##  rsvd                   1.0.3      2020-02-17 [1]
##  Rtsne                  0.15       2018-11-10 [2]
##  rvest                  0.3.6      2020-07-25 [1]
##  S4Vectors              0.24.2     2020-01-13 [2]
##  sandwich               2.5-1      2019-04-06 [2]
##  scales                 1.1.0      2019-11-18 [2]
##  scatterplot3d          0.3-41     2018-03-14 [1]
##  sctransform            0.2.1      2019-12-17 [2]
##  SDMTools               1.1-221.2  2019-11-30 [2]
##  sessioninfo            1.1.1      2018-11-05 [2]
##  Seurat               * 3.1.2      2019-12-12 [2]
##  SingleCellExperiment   1.8.0      2019-10-29 [2]
##  smoother               1.1        2015-04-16 [1]
##  sn                     1.5-4      2019-05-14 [2]
##  sp                     1.4-2      2020-05-20 [1]
##  stringi                1.5.3      2020-09-09 [1]
##  stringr              * 1.4.0      2019-02-10 [1]
##  SummarizedExperiment   1.16.1     2019-12-19 [2]
##  survival               3.1-8      2019-12-03 [3]
##  testthat               2.3.2      2020-03-02 [1]
##  TFisher                0.2.0      2018-03-21 [2]
##  TH.data                1.0-10     2019-01-21 [2]
##  tibble               * 3.0.4      2020-10-12 [1]
##  tidyr                * 1.1.2      2020-08-27 [1]
##  tidyselect             1.1.0      2020-05-11 [1]
##  tidyverse            * 1.3.0      2019-11-21 [2]
##  tsne                   0.1-3      2016-07-15 [2]
##  TTR                    0.23-6     2019-12-15 [1]
##  usethis                1.5.1      2019-07-04 [2]
##  uwot                   0.1.5      2019-12-04 [2]
##  vcd                    1.4-7      2020-04-02 [1]
##  vctrs                  0.3.5      2020-11-17 [1]
##  VIM                    6.0.0      2020-05-08 [1]
##  viridis              * 0.5.1      2018-03-29 [2]
##  viridisLite          * 0.3.0      2018-02-01 [2]
##  withr                  2.3.0      2020-09-22 [1]
##  xfun                   0.12       2020-01-13 [2]
##  xml2                   1.3.2      2020-04-23 [1]
##  xts                    0.12-0     2020-01-19 [1]
##  XVector                0.26.0     2019-10-29 [2]
##  yaml                   2.2.1      2020-02-01 [1]
##  zip                    2.0.4      2019-09-01 [1]
##  zlibbioc               1.32.0     2019-10-29 [2]
##  zoo                    1.8-7      2020-01-10 [2]
##  source                                 
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  Bioconductor                           
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Github (clauswilke/colorblindr@1ac3d4d)
##  R-Forge (R 3.6.2)                      
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  Bioconductor                           
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Github (saezlab/progeny@94bea1c)       
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.3)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
## 
## [1] /home/uhlitzf/R/lib
## [2] /usr/local/lib/R/site-library
## [3] /usr/local/lib/R/library